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ABSTRACT 

Acoustic  ray  paths  in  the  ocean  are  known  to  exhibit  sig- 
nificant horizontal  deflections  after  repeated  reflections 
from  the  bottom.   The  effect  may  be  quantitatively   and  quali- 
tatively observed  through  a  ray  trace  model  which  permits  a 
change  in  direction  of  the  vertical  plane  of  propagation  as  a 
function  of  bottom  slope  and  grazing  angle.   ECTRACE  is  a 
family  of  computer  programs  which  traces  a  bundle  of  rays  in 
three  dimensions  and  utilizes  bottom  depth  values  as  a  portion 
of  its  input.   Included  are  related  programs  which  develop  sea 
bed  models  from  digitized  bathymetry  data  or  synthetic  bathymetry 
functions. 
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I.   INTRODUCTION 

ECTRACE  is  a  FORTRAN  computer  program  designed  for  the 
investigation  of  bathymetric  effects  on  horizontal  acoustic 
ray  tracing.   Its  supporting  computer  programs  include  routines 
for  constructing  a  grid  model  of  the  sea  bed  from  contour  data, 
ray  tracing  routines,  and  associated  printed  and  plotted  output 
The  programs  were  written  for  IBM  360/370  series  systems  with 
a  Versatec  plotter  and  related  software.   The  significant  fea- 
ture of  the  tracing  algorithm  is  its  ability  to  shift  the  direc- 
tion of  the  vertical  plane  of  propagation  after  a  bottom 
reflection.   A  sea  bed  modeled  by  triangular  facets  fitted 
through  three  points  of  a  grid  cell  makes  this  feature  possible 
with  a  minimum  effect  on  computation  time. 

Several  models  for  tracing  acoustic  rays  in  the  ocean  are 
used  in  the  Navy  today.   These  models  are  based  on  a  simplified 
environment  because  of  practical  limitations  on  data  and  com- 
puter run  time,  and  they  adequately  trace  a  bundle  of  rays  in 
a  single  vertical  plane.   Some  enable  intensity  calculations 
at  a  given  distance  from  the  source.   TRIMAIN  [1],  for  example, 
is  a  model  which  features  range  as  well  as  depth  dependence 
of  sound  speed  and  accounts  for  some  interaction  with  irregu- 
lar bathymetry.   These  models  have  a  large  number  of  applica- 
tions in  deep  ocean  cases  involving  ducts,  channels,  and 
convergence  zones,  and  bottom  reflections  where  the  plane  of 
the  ray  path  does  not  exhibit  horizontal  deflection. 
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Actual  cases  of  long  range  propagation  over  irregular  bot- 
toms show  any  single  ray  path  rarely  remains  within  a  single 
vertical  plane.   Any  application  of  ray  theory  to  propagation 
involving  bottom  interaction  must  account  for  a  horizontal 
deflection  caused  by  ray  heading  changes  after  repeated  reflec- 
tion from  a  tilted  bottom  plane.   This  effect  is  most  pro- 
nounced for  long  range  paths  over  smoothly  sloping  surfaces, 
but  any  valid  ray  approximation  can  be  seen  to  undergo  signif- 
cant  heading  changes  after  only  a  few  reflections  from  an 
irregular  or  undulating  bottom. 

The  behavior  of  sound  propagation  over  troughs,  ridges, 
wedges,  and  seamounts  has  been  studied  analytically  by  Har- 
rison [2]  using  normal  mode  theory  and  ray  invariants.   Since 
actual  ocean  bottom  topography  over  large  area  defies  prac- 
tical analytical  description,  except  in  the  stochastic  sense, 
DeWitt  [3]  developed  numerical  techniques  for  three-dimensional 
ray  tracing  in  which  the  bottom  consists  of  triangular  facets 
generated  numerically  from  actual  or  interpolated  bathymetry 
data.   This  approximation  scheme  greatly  simplified  the  prob- 
lem of  calculating  the  three-dimensional  ray  path  parameters 
after  a  bottom  reflection  and  provides  a  means  of  constructing 
reasonably  accurate  ray  traces  (for  a  verification,  see  Appen- 
dix F)  . 

The  ECTRACE  modeling  program  is  an  adaptation  of  the 
Dewitt  technique  for  use  on  an  IBM  360/370  series  computer. 
The  initial  programs  may  be  used  to  generate  synthetic  bot- 
tom topographical  grids  or  create  approximations  to  real 
ocean  bottom  grids  from  a  bathymetric  data  source  file.   The 
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topographical  features  of  the  generated  grid  can  be  depicted 
by  a  contour  plot  and  a  perspective  surface  plot  for  compari- 
son with  the  original  contours.   The  primary  program  is  the 
ray  tracing  program  which  accepts  the  generated  bottom  grid 
and  performs  stepwise  ray  position  calculations  from  user- 
selected  initial  heading  and  elevation  angles.   The  output 
is  the  printed  history  of  each  ray,  and  a  two-dimensional  plot 
of  the  ray  paths  projected  on  a  vertical  plane  intersecting 
the  bottom  profile  along  the  axis  of  the  ray  fan.   Addition- 
ally, a  horizontal  plane  projection  is  plotted  which  may  be 
overlaid  on  the  appropriate  contour  plot  to  correlate  the 
horizontal  curvature  effect  with  grid  bathymetry.   Punched 
card  output  from  several  ray  tracing  runs  may  be  combined  in 
one  independent  routine  for  a  comprehensive  visual  inspection 
of  acoustic  shadow  zones  in  the  horizontal  plane. 

It  is  very  important  to  note  that  ECTRACE  permits  refrac- 
tion caused  by  sound  velocity  gradients  in  the  vertical  plane 
only  and  that  the  sound  velocity  in  the  horizontal  plane  is 
assumed  constant  within  a  well  defined  water  mass.   All  hori- 
zontal bending  effects  revealed  by  ECTRACE  are  the  result  of 
bottom  reflection  only. 

With  the  ability  to  convert  actual  bathymetry  data  to  a 
matrix  of  depth  values,  discrete  faceted  approximations  of  a 
selected  oceanic  region  can  be  used  for  generating  ray  traces 
of  operational  significance.   In  addition,  qualitative  studies 
of  horizontal  ray  deflections  over  idealized  bottom  configura- 
tions are  enhanced  by  the  computer-generated  plots  of  ray  paths 
over  easily  modifiable  synthetic  bathymetry.   ECTRACE  and  the 
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augmenting  programs  are  envisioned  to  serve  as  a  basis  for 
future  development  of  three-dimensional  ray  modeling  and  inves- 
tigation efforts. 
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II.   DESCRIPTION  OF  ECTRACE  AND  PERIPHERAL  PROGRAMS 

This  chapter  describes  each  of  the  programs,  their  pri- 
mary subroutines,  and  their  interrelationship  in  producing  a 
ray  trace  which  may  be  used  to  investigate  the  bathymetric 
effects  on  horizontal  ray  curvature.   They  are  presented  in 
their  intended  order  of  use.   Flow  Diagram  (Fig.  1)  and  Table 
I  are  provided  to  assist  in  understanding  program  relationships 

A.   GENBOT 

This  program  converts  contour  data  into  a  matrix  of  depth 
values  for  storage  on  a  permanent  device.   It  also  makes  a 
two  dimensional  contour  plot  of  the  input  data  and  super- 
imposes reference  latitude  and  longitude  lines  (Fig.  2). 
The  contour  plot  can  be  placed  underneath  the  horizontal  ray 
trace  plot  of  ECTRACE  or  ECCOM  on  a  light  table  to  aid  in 
the  visual  investigation  of  the  ray  curvature  (see  Figs.  3 
and  4,  and  Figs.  5,  6,  and  7). 

The  bathymetry  contour  data  used  in  this  research  exists 
as  a  sequence  of  data  sets  representing  ten  by  ten  degree 
regions  of  the  North  Atlantic  Ocean  and  were  obtained  exclu- 
sively from  a  tape  compiled  in  a  joint  project  of  the  Naval 
Research  Laboratory  and  the  Mobil  Oil  Company.   In  their 
original  form  the  data  were  sufficient  only  for  producing  a 
computer  plot  of  region  contours,  thereby  reproducing  por- 
tions of  a  chart  titled  "Bathymetry  of  the  Norwegian-Green- 
land and  Western  Barents  Sea"  [4]  (see  Fig.  8).   Their  use 
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in  GENBOT  required  us  to  insert  the  depth  values  of  the  con- 
tour lines  and  transfer  the  results  to  a  separate  tape  (de- 
scribed in  Appendix  A) . 

The  following  subroutines  are  called  by  GENBOT: 

1)  CORNER  calculates  boundary  latitudes  and  longitudes  for 
screening  only  those  data  needed  for  interpolation. 

2)  GEOPLT  superimposes  reference  latitude  and  longitude 
lines  on  the  subregion  contours. 

3)  GEODST  calculates  the  forward  or  inverse  solutions  of 
the  geodetic  triangle  to  accomplish  the  above. 

4)  RAIN  1,  RAIN  2,  and  RAIN  3  [5]  were  obtained  from  the 
SSP3  program  library  at  NPS  and  are  used  to  perform 
interpolations  on  data  points  for  construction  of  the 
depth  matrix. 

B .  SYNGEN 

This  program  produces  a  synthetic  bathymetry  grid  for 
storage  on  a  permanent  device.   It  may  be  used  to  generate 
the  following  types  of  ideal  sea  bed  configurations:   wedge, 
trough,  ridge,  conical  seamount  (Fig.  9),  sinusoidal  undula- 
tion, or  parabolic  basin  (Fig.  10).   These  grids  can  be  used 
in  the  ECTRACE  program  as  a  depth  matrix  in  the  same  manner 
as  the  grids  of  actual  bathymetry  produced  by  GENBOT.   SYNGEN 
is  discussed  in  detail  in  Appendix  B. 

C.  G3DP 

Program  G3PD  produces  a  perspective  of  any  portion  of  the 
generated  depth  matrix  through  use  of  the  NPS  system  subroutine 
CONTUR  [4].   Figures  9,  10,  and  11  are  examples.   G3DP  uses 
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job  control  commands  to  link  with  the  GRDSCT  subroutine  of  the 
ECTRACE  load  module  to  extract  the  desired  area  (hereafter 
called  the  working  matrix)  to  be  plotted.   Appendix  C  describes 
the  application  of  G3DP  in  detail. 

D .  GRDCHK 

This  program  is  used  to  check  the  validity  of  the  depth 
matrix  generated  by  the  programs  GENBOT  and  SYNGEN.   It  produces 
a  contour  plot  and  a  vertical  plot  of  user-selected  rows  of  the 
generated  matrix.   By  comparing  consecutive  rows  the  existence 
and  location  of  anomalous  data  points  may  be  determined  (Figs. 
4  and  12) .   It  contains  the  previously  mentioned  subroutine 
GEODST,  CORNER,  and  a  version  of  GEOPLT  (called  GEOPLR)  to  draw 
reference  geographic  lines  for  comparison  of  the  two  contour 
plots  (GENBOT  vs  GRDCHK) .   Refer  to  Appendix  D  for  additional 
discussion  of  GRDCHK. 

E .  ECTRACE 

The  ECTRACE  program  is  a  package  of  subroutines  which 
exist  in  load  module  form  on  a  permanent  storage  device  in 
the  NPS  computer  system  user  library.   Its  use  requires  a 
user  supplied  calling  program  which  calls  subroutine  TRACER 
and  contains  the  necessary  job  control  language  (JCL)  to  link 
with  the  module.   These  procedures  are  explained  in  more  de- 
tail in  Chapter  IV  of  this  report. 

ECTRACE  traces  rays  in  three  dimensions.   The  sound  speed 
field  is  assumed  to  vary  piecewise-linearly  with  depth  only, 
yet  provisions  are  made  to  permit  simulation  of  up  the  three 
distinct  water  masses  separated  discontinuously  by  vertical 
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fronts.   The  bathymetry  structure,  generated  by  one  of  the 
two  grid-making  programs,  is  a  three  dimensional  surface 
modeled  by  discrete  triangular  facets  fitted  through  cells 
of  adjacent  depth  values.   Rays  are  traced  as  piecewise  arc 
segments  each  with  a  radius  of  curvature  dependent  on  the 
vertically-structured  sound  speed  gradient. 

When  the  rays  undergo  a  bottom  bounce,  specular  reflection 
is  assumed.   Bottom  loss  values  can  be  calculated  from  the 
subprogram  FBLOSS  (described  later) ,  a  standard  Navy  model 
such  as  the  FACT  bottom-loss  routine,  or  a  user  supplied 
routine.   The  trace  history  supplied  by  ECTRACE  includes  the 
location  and  characteristics  of  each  ray  reversal  (surface 
reflections,  bottom  bounces,  and  refractive  turning  points). 
The  printed  data  lists  the  new  ray  parameters  determined 
after  the  reversal  had  occurred.   Chapter  V  explains  the 
printed  output  in  more  detail. 

The  ECTRACE  plot  product  includes  a  plane  view  of  the 
ray  paths  projected  on  a  vertical  plane,  a  plot  of  sound 
speed  profiles,  and  a  horizontal  (plane)  view  of  the  ray 
paths  (Figs.  5  and  13).   The  horizontal  view  can  be  overlaid 
on  the  appropriate  contour  plot  generated  by  GENBOT  or  GRDCHK, 
as  in  Figs.  12  and  14. 

The  following  subroutines  are  used  by  ECTRACE: 
1)   Tracer  is  the  primary  tracing  subroutine  which  calls 
all  other  subroutines  of  the  ECTRACE  program  and  is 
capable  of  processing  a  bundle  of  rays  separated  into 
ray  fans.   A  ray  fan  is  defined  here  as  a  number  of 
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rays  (user  determined)  each  with  a  different  initial 
heading  measured  clockwise  from  grid  north,  but  all 
with  the  same  initial  elevation  angle.   A  positive 
elevation  angle  describes  a  ray  which  is  pointed  toward 
increasing  depth. 

2)  GRDSCT  reads  into  core  a  fixed-dimension  section  of  a 
depth  matrix  file  of  arbitrary  size.   The  ray  trace  is 
restricted  to  the  boundaries  of  the  core-loaded  working 
matrix.   A  ray  which  departs  the  working  matrix  boundaries 
is  terminated  and  program  control  is  passed  to  the  next 
ray's  trace. 

3)  IDTSUB  identifies  the  three  matrix  subscript  pairs  whose 
depth  values  define  a  plane  triangular  bottom  facet 
underneath  the  ray  segment  head. 

4)  DEEP  calculates  the  depth  at  any  horizontal  position  by 
solving  for  the  z-coordinate  on  a  plane  fitted  through 
the  three  depth  matrix  values. 

5)  CONTAC  calculates  the  three-dimensional  coordinates  of 
the  intersection  of  the  ray  path  and  triangular  bottom 
facet. 

6)  BOANG  calculates  the  grazing  angle,  elevation  angle,  and 
new  heading  of  a  bottom  reflected  ray. 

7)  FBLOSS  returns  bottom  loss  values  in  dB  as  a  function 
of  grazing  angle  for  a  bottom  reflected  ray.   These 
values  are  based  on  NRL  geophysical  survey  data  of  the 
Greenland-Norwegian  Sea  and  Iceland  environment  collected 
in  the  mid-1970 's  and  have  been  widely  averaged  for 
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convenience.  The  user  has  the  option  of  substituting  a 
more  specialized  function  if  bottom  loss  values  are 
critical . 

8)  IDPROF  determines  the  water  mass,  and  hence  the  sound 
speed  profile,  affecting  the  ray  with  each  segment 
iteration. 

9)  NUPROF  initializes  calculations  for  new  water  mass  param- 
eters once  the  ray  has  passed  the  boundary. 

10)  CHNLIM  determines  in  advance  the  refractive  sound  channel 
limits  for  each  ray  from  its  invariant.   This  is  used  to 
test  a  ray  at  its  turn-around  point  and  to  reveal  its 
type  (surface  reflected,  ref racted-surface  reflected, 
bottom  reflected-surface  reflected)  in  the  printed 
trace  history. 

11)  ANGPRT  is  used  to  print  a  summary  of  selected  bottom 
contact  in  parameters  for  each  ray  fan. 

12)  The  following  subroutines  each  have  only  one  calling 
statement  in  TRACER.   The  subroutines  they  in  turn  invoke 
are  dependent  upon  the  computer  graphic  system  installed 
and  would  require  internal  modification  if  exported  to 

a  system  without  Versatec  software.   For  additional 
information  NPS  users  should  consult  Ref.  7. 

a)  BGNPLT  initializes  the  plotter  and  draws  all  borders, 
titles,  and  axes. 

b)  BDTPLT  draws  a  vertical  profile  of  the  sea  bed  cen- 
tered along  the  mean  heading  of  the  ray  bundle,  as 
shown  in  the  upper  plot  of  Figure  13. 
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c)  RNGPLT  traces  all  ray  segment  projections  on  the 
vertical  plot. 

d)  T2DPLT  plots  the  horizontal  track  of  each  ray.  Each 
ray  fan  uses  one  of  thirteen  symbols  to  indicate  the 
bottom  bounce  positions  of  each  ray.  Figure  5  shows 
some  of  the  available  symbols. 

e)  SSPPLT  draws  the  sound  speed  profiles. 

f)  ENDPLT  draws  the  horizontal  plot  legend  and  termi- 
nates all  plotting. 

F .   ECCOM 

Program  ECCOM  uses  the  optionally  produced  punched  card 
output  of  one  or  more  ECTRACE  jobs  to  draw  a  composite  hori- 
zontal plot  of  several  ray  bundles,  enabling  a  comprehensive 
pictorial  study  of  shadow  zones  and  acoustic  convergence  in 
the  horizontal  plane.   The  composite  horizontal  plot  is  de- 
signed to  portray  the  rays  which  emanate  from  the  same  source 
or  from  multiple  sources  along  a  desired  track.   The  ray  num- 
bers assigned  and  punched  for  each  ray  on  each  ECTRACE  run 
assist  the  user  in  singling  out  rays  of  interest  (or  non- 
interest)  before  combination.   Figures  7  and  18  are  examples 
of  this  combining  technique. 

Amplifying  remarks  on  these  programs  and  subroutines  can 
be  found  in  the  comments  in  the  computer  FORTRAN  listings. 
For  additional  description  of  plots  and  computer  printed  out- 
puts, refer  to  the  appropriate  appendices. 
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III.   PROGRAM  LOGIC  AND  THEORY 

ECTRACE  traces  the  three-dimensional  paths  of  a  bundle  of 
rays,  one  ray  at  a  time.   The  bundle  is  divided  into  fans  of 
rays  of  the  same  initial  elevation  angle.   A  separate  ray  fan 
is  traced  for  each  fixed  increment  of  initial  elevation  angle 
between  the  limits  ELST  and  ELEND  (input  constants) .   A  posi- 
tive elevation  angle  describes  a  ray  pointed  downward.   The 
initial  ray  headings  in  each  ran  fan  extend  from  KDST  to  HDEND 
(input  constants) ,  measured  in  degrees  clockwise  from  grid 
north. 

A.   BOTTOM  AND  SOUND  FIELD  MODELING 

Figure  III.l  shows  the  projection  of  the  bottom  surface 
onto  the  horizontal  x-y  plane. 
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FIGURE  III.l   Bottom  Surface  Projection 
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The  projection  of  each  triangular  facet  is  a  right  isosce- 
les triangle.   The  entire  surface  can  be  defined  by  specifying 
the  depth  values  at  the  vertices  of  the  projection  triangles. 
ECTRACE  stores  these  depth  values  (km)  in  the  two-dimensional 
working  array  ZB.   The  parameter  DHN  (km)  represents  the  spac- 
ing between  the  vertices  in  both  the  X  (East)  and  Y  (North) 
directions.   Thus,  the  value  of  the  matrix  point  ZB(I,J)  is 
the  depth  at  the  vertex  point  with  X  coordinate  (I  -  1)  * 
DHN  and  Y  coordinate  (J  -  1)  *  DHN,  relative  to  the  ZB  origin. 

The  matrix  values  for  ZB  are  read  from  a  source  file  stored 
on  a  permanent  storage  device  such  as  a  disk  pack.   A  matrix  of 
dimensions  151  by  151  with  spacing  DHN  =  km  has  been  used  for 
trial  runs. 

1.   Bottom  Contact  Point 

The  subroutine  CONTAC  locates  the  coordinates  and  calcu- 
lates the  depth  of  the  bottom  bounce  point.   The  parameters 
describing  locations  of  the  ray  segment  head  and  tail  are  among 
those  passed  to  the  main  program.   In  this  subroutine,  the  ray 
segment  is  treated  as  unrefracted  (a  straight  line) ,  a  reasonable 
approximation  for  the  small  distances  involved  in  determining 
the  contact  point. 

First  the  horizontal  positions  of  the  ray  segment  head  (x,y) 
called  CEE  and  CNN  in  CONTAC,  and  (xt,y2)  and  called  CRE  and 
CRN  in  CONTAC,  are  checked  to  see  if  they  both  lie  in  the  same 
projection  triangle.   If  not,  an  interative  routine  searches 
for  two  points  on  the  segment  such  that  the  horizontal  coordi- 
nates of  both  points  lie  inside  the  same  triangular  cell  and 
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define  the  portion  of  the  segment  that  penetrates  the  bottom 
facet.   Then  the  depth  of  the  ray  head  z,  and  the  depth  of  the 
tail  zl    are  computed. 

Figure  III.  2.  illustrates  a  profile  of  a  cell  bottom  with 
a  ray  intersecting  it  in  the  propagation  plane.   The  solution 
of  the  triad  (x  ,y  ,z    )  identifying  the  point  of  bottom  contact 
proceeds  as  follows: 

Using  subroutine  DEEP  (described  in  Section  3)  compute  H   and 
H2,  the  bottom  depth  at  (x  ,y  )  and  (x,/y2)/  respectively.   Then 
compute : 

Az  =  z   -  z 

1  2 

where  (x  ,y  ,z    )  are  contact  point  parameters. 

C     w     w 

Hz   -  Hz 

Z    =    1   2 2  1 

C   Az  +  Ha  -  H2 


If      AZ  ¥■    0,  M  =   fc Z_j_ 

AZ 


If      AZ  =  0,  M  =   2i  "  Hi 

Hi  ~  H2 


Then   x  =  x  +  M  r  sin  <j>     and 
c 


y   =  y   +  M  r  cos  (J> 
where  <J>  is  the  ray  heading  (angle  between  the  y-z  plane  and  the 
propagation  plane) . 

2.   Projection  Triangle  Vertex  Subscripts 

The  subroutine  IDTSUB  calculates  the  I  subscripts,  IRT, 
ILT,  IVR  and  the  J  subscripts  JBT,  JTP ,  JVR  for  the  vertices  of 
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FIGURE  III.  2.   Profile  of  the  bottom  cell  in  the  propagation 

plane,  with  a  ray  contact. 


of  the  projection  triangle  which  contains  the  ray's  horizontal 
coordinates  x  and  y.   The  subscript  parameters  are  indicated 
in  Figure  III. 3.   The  following  suffixes  identify  the  relative 
locations  of  the  vertices  of  the  projection  triangle: 

RT  -  right 

LT  -  left 

TP  -  top 

BT  -  bottom 


iLTjOTf  * 


tuTjTaT 


KIT,  Iff* 


•*Xfc,T,TBf 


Xv/e,l\/fcf 


xaT,  Trr 


1  LT j  T&T 


TLTj  T6T" 


X&T,  TT? 


*  IVfc,  J"v& 


FIGURE  III. 3.   Projection  Triangle  Vertices 
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3.   Depth  Calculation 

Subroutine  DEEP  calculates  the  depth  in  kilometers 
corresponding  to  the  horizontal  coordinates  x  and  y.   Auto- 
matically, IDTSUB  is  called  to  identify  the  projection  tri- 
angle subscripts.   The  depth  z   is  calculated  by  solving  the 
equation  of  a  plane  defined  by  the  depths  at  the  vertex  points 
(D^  D  ,  and  D  ).   Figure  III.  4  illustrates  the  calculations 
made  based  on  the  orientation  of  the  projection  triangle  out- 
lining the  bottom  facet  plane.   The  quantities  x   and  y   repre- 
sent the  horizontal  position  (in  units  of  DHN)  of  the  depth 
calculation  point  relative  to  the  horizontal  position  of  the 
depth  calculation  point  relative  to  the  horizontal  position 
of  the  depth  value  D1     (the  lower  left  vertex  of  the  projec- 
tion triangle)  as  calculated  in  Secion  1.   The  depth  values 

D   and  D  '  are  normalized  by  the  length  of  the  equilateral 
x      y  * 

side  1,  in  units  of  DHN.   Finally,  the  depth  is  solved  by 


z   =xD   +yD   +  D, 
c     ex     c  y     i 


B.   RAY  LOCATION  AND  DIRECTION  PARAMETERS 

A  three  dimensional  left  hand  coordinate  system  is  used 
(Fig.  III. 5.),  with  the  positive  z  axis  pointing  in  the  direc- 
tion of  increasing  depth  and  the  sea  surface  level  lying  in 
the  x-y  plane.   As  the  ray  is  traced  through  three  dimensional 
space,  the  program  variables  which  are  used  to  locate  the  tail 
of  the  ray  vector  are: 

DEPl  -  depth  coordinate 

CRE   -  x  coordinate  in  km 
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FIGURE  III. 4.   Depth  calculation  on  bottom  facet 

CRN  -  y  coordinate  in  km 
CRE  and  CRN  are  components  of  the  horizontal  range  increment 
DR  which  forms  the  propagation  plane  heading  angle  PHI  with 
the  y-axis. 
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FIGURE  III. 5.   Ray  Location  Parameters 
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Given  that  the  ray  vector  originates  at  the  point  with 
the  above  coordinates,  the  direction  of  the  ray  vector  is 
specified  by  the  following  angles: 

TH1  -  ray  elevation  angle  in  the  plane  of  propagation 
measured  from  the  ray  vector's  projection  on  the  xy 
plane  to  the  ray  vector.   A  ray  pointing  toward  the 
ocean  bottom  will  have  a  positive  elevation  angle. 
Theta  takes  on  values  between  -tt/2  and  +t/2  (output 
values  are  in  degrees) . 

PHI  -  ray  (propagation  plane)  heading  angle  between  the 
positive  y-axis  and  the  projection  of  the  ray  vector  on 
the  xy  plane,  measured  clockwise  from  the  positive  y-axis 
(grid  North).   PHI  takes  on  values  between  -it  and  It. 
Initial  values  for  THl  and  PHI  are  calculated  from  user  input 
values  (in  degrees) ,  which  are  converted  to  radians  for  compu- 
tational purposes.   The  unit  vectors  in  the  horizontal  plane 
are  given  as  UVX  =  sin  (PHI)  and  UVY  =  cos  (PHI) . 

C.   RAY  PATH  CALCULATIONS 
1.   Path  Length 

ECTRACE  confines  refraction  to  the  vertical  plane  alone. 
It  is  recognized  that  while  horizontal  sound  speed  gradients 
exist,  their  effects  are  sufficiently  small  within  a  well  defined 
water  mass  to  be  neglected  compared  to  the  bathymetric  effect. 

However,  while  the  horizontal  gradients  are  likely  to 
be  slight,  important  discontinuities  in  water  characteristics 
may  occur  horizontally  across  oceanic  fronts.   For  this  reason, 
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ECTRACE  has  been  designed  to  simulate  up  to  three  distinct 
water  masses  and  two  associated  frontal  boundaries. 

Within  a  water  mass  having  uniform  characteristics 
along  fixed  depths  a  ray  path  conveniently  defines  a  vertical 
plane  within  which  Snell's  Law  can  be  written  as 

V  =  c(z)/cosO  (z)  (1) 

where  V  is  treated  as  invariant  within  the  water  mass  (termed 
the  vertex  velocity  in  some  texts  since  it  is  the  speed  at  the 
depth  at  which  the  ray  vector  becomes  horizontal) .   The  term 
c(z)  is  the  sound  speed  at  depth  z  and  6  (z)  is  the  elevation 
angle  at  depth  z. 

Equation  (1)  uses  the  simplification  that  the  speed  of 
propagation  and  the  elevation  angle  are  functions  of  depth  only 
It  is  generally  possible  to  define  a  stratified  medium  in  which 
each  layer's  gradient  is  constant.   Although  such  a  model  re- 
quires a  large  number  of  data  points  to  adequately  describe  a 
water  column  of  complicated  sound  speed  structure,  the  piece- 
wise  linear  gradient  approximation  greatly  simplifies  the  ray 
path  descriptions. 

Taking  advantage  of  the  approximations,  we  use  a  simple 
equation  for  the  sound  speed  gradient,  within  a  layer, 

g  =  dc/dz  (2) 

Integrating  yields 

z  =  c/g  (3) 
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where  the  constant  of  integration  is  avoided  by  arbitrarily 
specifying  the  origin  of  the  coordinate  axes  at  the  depth  where 
c  =  0.   By  applying  Snell's  Law  we  obtain 

z  =-—  cos6  (4) 

and 

V     „  „ 
dz  = sinGdS  (5) 

Specifying  r  as  the  horizontal  distance  axis  of  the  vertical 
plane  and  applying  analytical  geometry  there  results 

dz  =  tanGdr  (5) 

dr  = cos6de  (7) 

and 


r  -  -  JL.  sine  (8) 

when  measuring  r  from  the  vertical  (z-axis)  to  the  point  at 
which  9=0  (horizontal  projection) . 

Squaring  (4)  and  (8)  and  summing  yields 

r2  -  z2  =  C-^)2  (9) 

This  is  the  equation  of  a  circle  of  radius  (V/g)  and  whose 
center  is  at  the  origin  of  the  coordinate  system  just  described 
The  equation 

R \-  (10) 

describes  the  radius  of  curvature  of  a  ray  path  in  a  constant 
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sound  speed  gradient.   The  minus  sign  has  been  chosen  to  allow 
3  to  be  positive  when  the  refracted  ray  is  increasing  its  ele- 

vation  angle. 

The  coordinate  system  and  computer-approximated  ray 

path  segments  are  shown  in  Fig.  III. 6. 
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FIGURE  III.  6.  Ray  Path  Geometry 
||         it  can  be  shown  from  Fig.  III. 6.  that  the  ray  segment 
path  length  between  P  and  Q  is 

A£  =  pA6 


til) 


where 


A0  = 


(12) 


ion  angle  of  the  ray  is  presumed  known, 


Since  the  initial  elevat 

dues  of  6  are  found  iteratively: 


successive  va. 


Ail 


(13) 
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L  geometry  also  allows  for  solution  of  the  depth  change, 
iz  .  z   .  z   =  p(cos8   -  oos62)    ,      C14) 

2       1 

L  likewise  the  increase  in  horizontal  range, 

Ar  .  ^  .  tj  =  p(sin62  -  sin9i)   .      (15) 

L  that  in  the  «g.  XXX...  «■  —  —  -  P"**"' 

0  the  sign  of  P  ■»*  agree  with  that  of  .a.   Sguatio n  (1 

hich  gives  p  the  opposite  sign  of  the  gradient,  satrsfre  th  s 

equirement.   Only  Az  is  allowed  to  heco.e  negative,  as  when 

ray  decreases  its  depth. 

2   Travel  Time 

'  cvia      III  7  )  the  sound  speed  relation- 

From  the  diagram  (Fig.  Ill- /.J 


ship  becomes 


c  = 


df 
dt 


(16) 


(17) 


>  r 


FIGUR^  III. 7.   Ray  path  and  s 


ound  speed  relationship 
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The   ray   invariant 


.-i-22S£  U8) 


and  the  trigonometric  identify  sin0  =   \jl  -  cos2  6    combine 
and  give  dt  as  a  function  of  depth 


dt  =  dz      (19) 

c(t)  vli  "  Uc(z)]2 


f~- 


The  time-of-f light  of  the  acoustic  ray  along  the  path  is  found 
by  integrating  from  the  initial  point  to  the  final  point: 


1  dz 

t^  -  t.  =   / 


f    ±  «±   c(z)  f  -  [ac(z)]2     (20) 

Before  integrating  it  is  assumed  that  the  gradient  is  constant 
along  the  segment.   More  generally,  ECTRACE  assumes  a  constant 
gradient  in  the  layer  bounded  by  z   and  z  . 
Then 


c(z)  =  c(z  )  +  g(z   -  z  )  (21) 

1  2       1 

where  z  <  z  <,    z 

1  2 

The  integral  may  now  be  evaluated  by  introducing  a  new  variable 
W  =  £iil  such  that 

g 

w  =  z  _  z   +  ^i-  (22) 

1    g 

Then  dW  =  dz  and  C (W)  =  gW.   Therefore, 

Wf        dw 


t,  -  t.  =   / 

r       l        VJ;     _T.T   .h      ^2„2TT2 


^i   gW  \|l  -  a2g2w: 
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1      W,[l  +  P  -  a2g2Wc2]        (23) 

-  £n   r    t » . } 

9      w. [1  +  I  1  "  azgzWf2] 


W   CI  +  sine  . ) 
„  1  to  [  !lJ i-  ]  C24) 

9      W   (1  +  sine,;) 
1  r 

Values  for  W  are  restricted  to  being  in  the  s»e  layer  and  on 
the  same  side  of  a  turning  point  in  the  ray  path. 

D.   RAY  TRACING  PROCEDURE 

*     i  i.  r-   ni,i0n  initial  coordi- 
ECTRACE  begins  tracing  a  ray  from  its  given 

nates  CRE ,  CRN,  and  DEPl.   First  the  water  mass  is  identified 

to  determine  the  local  gradient  and  therefore  the  ray  invariant 

ThP  rav  oath  is  constructed  in  the 
and  radius  of  curvature.   The  ray  patn 

._  i.rrpment  DTH  calculated  from  a 
form  of  arc  segments  using  an  increment  u 

chosen  arc  length  «..      Then  the  new  6  value  is  calculated  and 
nations  C14)  and  (15,  are  used  to  determine  the  ray  segment's 
aepth  increment  and  horizontal  range  increment.   The  horizontal 
position  is  updated  using  the  horizontal  unit  direction  vectors 

-,  *     ntiT   rhp  initial  heading  of  the 
(UVX  and  UVY)  calculated  from  PHI,  the  imt 

ray's  plane  of  propagation. 

a»*e  3-D  coordinates  of  the  ray  head  may  now 
The  interim  update  3-D  cooiui 

be  stated  in  FORTRAN  as 

(25) 
CEE  =  CRE  +  DR  *  UVX 

(26) 
CNN  =  CRN  +  DR  *  UVY 

(27) 
DEP2  =  DEPl  and  DZ 

«.*■  the  horizontal  range  increment  Ar  and  DZ  is 
where  DR  represents  tne  non 

the  depth  incrementAz. 
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E.   RAY  CONSTRAINT  CONDITIONS 

Each  time  the  ray  path  is  incremented,  certain  constraint 
conditions  must  be  checked.  These  constraints  are  bottom  re- 
flections and  the  sound  channel  limits. 

1 .   Bottom  Reflections 

After  the  interim  horizontal  coordinates  CEE  and  CNN 
are  calculated,  they  are  passed  to  DEEP  to  find  the  bottom 
depth  WH2  at  the  ray  head.   If  DEP2  is  greater  than  or  equal  to 
WH2 ,  a  bottom  bounce  has  occurred.   In  this  case  the  values  of 
DEP2,  CEE,  and  CNN  are  adjusted  by  linear  interpolation  and  a 
point  of  bottom  contact  is  established.   Subroutine  BOANG  is 
then  called  to  calculate  new  values  of  TH1,  PHI,  UVY,  and  the 
ray  grazing  angle  GRAZ . 


bottbrA  -foce? 


Plan€  of 


FIGURE  III. 8.   Bottom  Bounce  Reflection  Angles. 

The  basic  calculations  of  BOANG  involve  the  vectors: 

v  -  unit  ray  vector   before  bounce 

n  -  unit  vector  normal  to  triangle  facet  at  the  bounce  point 

v'-  unit  ray  vector  after  bounce 
The  components  of  v  are  calculated  from  TH1  and  PHI.   The  com- 
ponents of  n   are  calculated  from  the  equation  of  the  plane  of 
the  bottom  facet.   The  new  vector  v  is  then  calculated  from 
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the  vector  equation, 

v'  =  v  -  2  (v  •  n)  n    .  (28) 

New  values  of  TH1,  PHI,  UVX  and  UVY  are  then  calculated  from 
the  compnents  of  v'.   Finally,  the  grazing  angle  GRAZ  comes 
from 

sin  (GRAZ)  =  -Cv  •  n)   .  (29) 

Reflection  from  a  sloping  bottom  plane  facet,  as  op- 
posed to  a  surface  reflection,  forces  a  reevaluation  of  the 
Snell ■ s  Law  constant  V,  and  therefore  the  radius  of  curvature 
R.   The  invariant  V  remains  in  effect  until  the  next  bottom 
bounce  occurs  or  a  new  water  mass  is  entered.   Thus,  after 
summoning  BOANG,  the  ray  tracing  continues  with  the  new  ray 
parameters.   The  sound  channel  limits  CHS  and  CHD  described 
below  are  also  recalculated  at  this  point  since  they  also 
depend  on  V,  and  GRAZ  is  passed  to  function  FBLOSS  to  calcu- 
late a  bottom  loss  value. 

2 .   Sound  Channel  Limits 

The  Snell' s  Law  equation  (1)  says  that  the  ray  trace 

is  constrained  to  depths  where  the  speed  of  sound  does  not 

exceed  V.   A  given  ray  path  may  be  constrained  in  the  sense 

that  it  is  purely  refracted  (R) ,  refracted  or  reflected  (RSR 

or  RBR)  or  purely  reflected  (BRSR) .   The  determining  factors 

are  V  and  the  local  sound  speed  by  profile.   V  is  calculated 

from  the  sound  speed  at  the  initial  depth  and  the  initial 

elevation  angle,  and  is  recalculated  after  a  bottom  bounce  or 

upon  entering  a  new  water  mass. 
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CHS  and  CHD  are  the  refractive  sound  channel  limits 
at  the  shallow  and  deep  extents  respectively.   Since  the  sound 
speed  profile  is  approximated  by  piecewise  continuous  linear 
gradients,  the  water  column  is  treated  as  a  column  of  indexed 
layers,  each  characterized  by  upper  and  lower  (shallow  and 
deep)  depth  boundaries  and  a  constant  sound  speed  gradient. 
The  limiting  depths  CHS  and  CHD  are  found  by  scanning  the 
water  mass  layers  both  above  and  below  the  current  ray  loca- 
tion until  a  boundary  sound  speed  value  is  found  which  exceeds 
V.   The  limiting  depth  is  then  calculated  from 

z,  .   =   z  +  (V  -  c)/g  (30) 

lim  ' 

where  z  and  c  are  evaluated  at  the  top  of  the  layer  and  g  is 
the  gradient  of  the  limiting  layer. 

For  cases  where  g  =  0  or  where  V  is  not  exceeded  by 
any  sound  speed  in  the  local  profile,  CHS  and  CHD  are  set  to 
the  physical  limits  of  the  water  column. 

F.   RAY  PROPAGATION  TERMINATION  CONDITIONS 

Tests  are  made  during  each  ray  segment  iteration  for  events 
which  terminate  the  ray  trace.   These  events  are: 

•  Escape  from  grid  -  when  the  ray  departs  the  finite  x  and 
y  constraints  of  the  working  depth  array  (ZB) . 

•  Depth  boundaries  -  when  the  depth  detected  at  the  ray  is 
less  than  ten  wavelengths  or  greater  than  10  km. 

•  Bottom  loss  -  when  the  attenuation  due  to  bottom  reflec- 
tion exceeds  a  specified  input  amount  in  dB. 
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*  Error  at  bottom  contact  -  when  defective  program  code 
or  anomalous  depth  data  prevent  the  interpolative  rou- 
tine CONTAC  from  recovering  the  bottom  depth  within  a 
reasonable  number  (50)  of  iterative  attempts. 

*  Trapping  -  when  it  becomes  apparent  that  the  ray  will 
not  contact  the  bottom,  rendering  further  three-dimen- 
sional tracing  unnecessary. 

*  Bottom  bounce  gate  -  when  the  total  number  of  bottom 
reflections  exceeds  a  user  specified  amount. 

The  cause  of  ray  termination  is  specified  in  the  printed 
ray  history,  giving  the  ray  parameters  at  the  time  of 
termination . 

G.   IMPROVEMENTS  MADE  OVER  THE  ORIGINAL  NRL  VERSION 

The  ray  tracing  algorithm  originally  existed  as  a  program 
called  ABOUNC  at  the  Naval  Research  Laboratory.   ECTRACE  dif- 
fers significantly  from  ABOUNC  in  the  following  areas: 
1 .   Ray  Path  Development 

ABOUNC  calculates  ray  parameters  through  right  tri- 
angle solutions  and  Snell ■ s  Law,  but  develops  the  ray  path 
along  fixed  depth.   Since  this  prevents  rays  from  becoming 
horizontal,  a  reversal  is  forced  when  the  magnitude  of  the 
elevation  angle  goes  below  a  gate  value.   In  contrast,  ECTRACE 
develops  the  ray  path  along  arc  segments  and  radii  of  curva- 
ture, fixing  only  the  maximum  length  of  the  segments  while 
interrupting  the  increment  as  necessary  when  a  segment  reaches 
a  reversal  point  or  layer  boundary,  or  to  keep  the  increment 
within  a  maximum  amount. 

39 


2 .  Computer  Word  Storage 

ABOUNC  floating  point  variables  are  nearly  all  FORTRAN 
double  precision.   The  largest  single  burden  on  computer  core 
is  the  depth  matrix,  read  from  a  magnetic  tape  in  card  image 
format,  each  card  record  containing  eight  items  (fields)  ten 
bytes  long.   This  format  permitted  a  depth  range  of  ±9999  km 
with  resolution  of  ±1  cm.   By  reducing  the  depth  range  to  ±999 
km  and  resolution  to  ±10  cm,  two  unnecessary  bytes  are   removed 
from  each  field,  enabling  a  fit  of  ten  items  per  record  and  a 
permanent  storage  reduction  of  twenty  percent.   In  addition, 
ECTRACE  reads  the  matrix  into  a  single  precision  array  to  re- 
duce the  core  requirements  for  the  data  by  fifty  percent. 

3.  Grid  File  Manipulation 

ABOUNC  reads  the  contents  of  the  entire  depth  matrix 
data  file  into  array  v/hose  dimensions  must  be  exactly  equal  to 
those  of  the  matrix.   This  requirement  is  not  inconvenient  on 
a  large  computer  with  FORTRAN  dialect  which  permits  objects 
time  dimensioning.   ECTRACE  is  written  in  FORTRAN  G  for  an  IBM 
Series/360  Model  67  machine  whose  limitations  required  some 
fundamental  modifications  of  the  program  and  its  job  control. 

To  maintain  flexibility  and  reduce  run  time,  ECTRACE 
was  reduced  to  a  package  of  subprograms  which  have  been  pre- 
compiled and  stored  on  a  permanent  device  in  load  module  form. 
The  user  must  supply  a  calling  program  which  serves  primarily 
to  establish  the  dimensions  of  the  depth  array.   The  user's 
selection  of  array  dimensions  are  made  on  the  basis  of  core 
economy  and  area  of  interest  rather  than  the  dimensions  of  the 
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input  matrix.   ECTRACE,  through  its  internal  subroutine  GRDSCT, 
extracts  a  rectangular  working  subregion  from  any  part  of  the 
input  matrix.   Trial  runs  of  ECTRACE  used  a  working  array  of 
dimensions  151  by  151  representing  a  sector  of  150  km  on  a  side 
with  a  core  requirement  of  less  than  250K,  while  the  input  ma- 
trix was  36  9  by  4  43.   Without  GRDSCT,  ECTRACE  would  have  re- 
quired an  allocation  of  800K  bytes  using  the  same  input  file. 
4 .   Treatment  of  Small  Elevation  Angles 

As  stated  previously,  ABOUNC  forced  a  reversal  of  the 
ray  as  the  elevation  angle  of  the  leading  segment  approached 
a  gate  value  close  to  zero.   Since  the  ECTRACE  tracing  scheme 
permits  a  full  range  of  elevation  angles,  early  trial  runs 
revealed  an  inherent  error  phenomenon  as  0  approached  zero  in 
a  steep  gradient.   Specifically,  a  A6  on  the  order  of  3  de- 
grees caused  by  a  small  radius  of  curvature  or  along  arc  seg- 
ment traced  a  ray  to  its  vertex  depth  in  error  as  much  as  50 
meters  from  the  Snell's  Law  prediction.   Since  the  depth  change 
was  calculated  from  the  difference  in  elevation  angle  cosines, 
the  error  was  found  to  arise  from  the  machine  inability  to 
retain  difference  precision  for  any  pair  of  numbers  close  to 
integer  values  (round  off  error) .   Since  the  difference  rather 
than  the  actual  cosine  values  is  needed  for  this  calculation, 

Q  2        Q    4 

the  small  angle  series  approximation  cos     1  -  y  +  -j,  -  .  .  . 
was  used  to  formulate  the  equation 

cosO   -  cos0  =   h(Q2-Qz)    -  HS"    -Bk).  (31) 

1  2  2     1       24    2      1 

When  restricted  to  small  angles,  this  approximation  yields 
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greater  accuracy  than  the  straightforward  machine  calculation 
of  the  left  hand  side,  since  the  precision  retained  by  the 
registers  is  very  high. 

5 .  Travel  Time 

ABOUNC  calculates  the  travel  time  of  an  acoustic  sig- 
nal along  a  ray  segment  by  dividing  the  linearly-approximated 
path  length  by  the  mean  sound  speed.   ECTRACE  accounts  for 
path  curvature  and  depth  dependence  of  the  sound  speed  through 
integral  evaluation,  discussed  in  detail  in  section  C.2. 

6 .  Track  Events 

Layer  transitions,  reflections,  and  turnarounds  are 
events  which  may  require  a  change  in  ray  parameters.   ABOUNC 
tests  for  these  events  in  algorithmic  sequence,  makes  the 
necessary  changes  in  response  to  the  first  event  detected 
and  continues  the  trace  from  the  current  position  of  the  ray 
vector  head.   ECTRACE  retraces  a  ray  segment  to  the  event 
depth  and  tests  for  all  other  events  before  beginning  the  next 
segment  from  the  transition  point.   This  logic  eliminates  the 
possibility  of  a  long  ray  segment  crossing  more  than  one  event 
depth  with  only  one  event  detected. 
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IV.   USER  INSTRUCTIONS  (ECTRACE) 

A  single  step  ECTRACE  job  traces  a  bundle  of  rays  from  a 
single  source  in  the  form  of  ray  fans.   Each  ray  fan  is  a  bun- 
dle of  rays  of  common  initial  elevation  angle.   All  ray  fans 
are  bounded  between  to  initial  heading  values  measured  from 
grid  North  (grid  North  equals  true  North  at  the  grid  center) . 

During  the  following  discussion  the  reader  may  wish  to 
refer  to  Fortran  listing  in  Appendix  M.   The  user  instructions 
for  the  other  programs  (GENBOT,  G3DP,  SYNGEN,  GRDCHK,  and 
ECCOM)  are  explained  in  their  individual  appendices  (User  In- 
structions and  Output  Description) . 

A.   JOB  CONTROL  CARDS  (JCL) 

The  following  JCL  cards  and  parameters  are  required  for 
an  ECTRACE  run: 

•  JOB  card  -  The  CPU  time  parameter  should  allow  60  sec- 
conds  for  every  seven  rays  to  be  traced. 

•  EXEC  card  -  This  card  must  specify  the  FORTCLGW  proce- 
dure (Appendix  0) .   250K   bytes  of  core  is  required  for 
the  standard  22801  (151  by  151)  element  depth  matrix 
portion  plus  4K  for  every  additional  1000  elements. 
The  computer  systems  job  class  definitions  should  be 
considered  before  deciding  upon  the  size  of  the  working 
matrix  and  the  number  of  rays.   For  example, 

EXEC  FORTCLGW, REGION. GO=250K 
is  a  class  C  job  at  NPS  when  limited  to  five  minutes  of 
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time,  which  is  adequate  for  tracing  4  0  rays. 

•  Calling  program  -  Following  the 

//FORT. SYS IN  DD  * 
card,  the  calling  program  must  contain  cards 

COMMON/DIMS/KMX , KMY  and 

CALL  TRACER  (ZB)  . 
ZB  is  a  REAL  *  4  array  dimensioned  exactly  by  KMX,KMY,  and 
these  values  (KMX,KMY)  are  initialized  to  the  values 
of  the  DIMENSION  statement  using  assignment  statements, 
as  in  the  following  example. 

DIMENSION  ZB(151,151) 

COMMON  /DIMS/  KMX,KMY 

• 

KMX=151 

KMY=151 

CALL  TRACER (ZB) 

STOP 

END 
KMX  and  KMY  are  INTEGER  *  4  and  should  be  less  than  or 
equal  to  the  dimensions  of  the  data  matrix  from  which 
the  depth  values  will  be  taken. 

•  LINK  cards  -  These  cards  point  to  the  ECTRACE  load 
module.   They  must  include  all  members  of  the  module 
for  which  the  user  does  not  supply  a  substitute.   An 
example  which  includes  all  members  is: 

//LINK.USDD  UNIT=3  33  0,VOL=SER=DISKOl,DSN=127  0.ECTRACE, 
//   DISP=SHR,LABEL=( , , ,IN) 
//LINK.SYSIN  DD  * 
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INCLUDE  USDD (TRACER, GRDSCT , IDPROF ,NUPROF , CKNLIM) 
INCLUDE  USDD (BGNPLT , SSPPLT , BDTPLT , RNGPLT , T2DPLT , ENDPLT ) 
INCLUDE  USDD ( IDTSUB , DEEP , CONTAC , BOANG , ANGPRT , FBLOSS ) 
ENTRY  MAIN 

/* 

•  GO  cards  -  ECTRACE  uses  data  set  reference  number  one 
for  input  of  the  depth  matrix.   For  example,  to  point 

to  a  depth  matrix  residing  in  a  data  set  called  S9999.GRID 
on  DISK03,  the  following  DD  statement  is  used: 

//GO.FT01F001  DD  UNIT=330 , VOL=SER=DISK03 , DSN=S9999 . 
GRID, DISP=SHR, LABEL= ( , , ,IN) 
If  a  depth  matrix  from  an  outside  source  is  used,  it 
must  be  in  the  format  specified  in  the  Output  Grid  Data 
Format  section  of  Appendix  A.   If  a  punched  card  output 
is  not  desired  to  make  an  ECCOM  run,  include  the  card 

//GO.FT07F001  DD  DUMMY 
ahead  of  the  GO.FT01F001  DD  statement  shown  above. 

B.   DATA  CARDS 

Following  the 

//GO.SYSIN  DD  * 
card,  the  user  selected  data  must  be  provided  in  the  following 
format.   FORTRAN  field  descriptors  are  in  parentheses.   The 
FORTRAN  variable  names  are  listed  for  reference. 

•  Card  1  (20A4) 

OTIT  -  title  of  the  run 

•  Card  2  (958.2) 

XO  -  source  position  for  x  coordinate  (kilometers) . 
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YO  -  source  position  y  coordinate  (kilometers) . 
ZO  -  source  depth  (meters) .   A  source  depth  which 

exceeds  the  water  depth  will  be  decreased  to  one 
meter  above  the  bottom.  Source  coordinates  in  a 
nonpropagating  location  will  also  be  adjusted. 

XREC  -  receiver  position  x  coordinate  (kilometers) . 

YREC  -  receiver  position  y  coordinate  (kilometers) . 

ZREC  -  receiver  depth  (meters) .   The  receiver  may  be 

positioned  at  the  source  if  calcuations  are  not 
desired.   Excessive  depth  is  automatically 
corrected. 

FREQ  -  frequency  (Hz) .   A  ray  is  terminated  when  the 
water  depth  becomes  less  than  ten  wavelengths 
during  propagation. 
FDL  -  (optional)  segment  length  (meters) .   The  parameter 
affects  the  smoothness  of  the  vertical  plot. 
Smaller  values  increase  smoothness  and  run  time. 
Large  values  create  the  risk  of  missing  a  bottom 
peak.   The  default  value  is  one  wavelength. 
•  Card  3  (6F8.2) 

XMIN  -  initial  x  value  of  horizontal  plot  (kilometers) . 

YMIN  -  initial  y  value  of  horizontal  plot  (kilometers) . 

The  above  values  are  used  for  drawing  the  axis  of  the 

horizontal  plot  and  for  determining  the  section  of  the 

input  matrix  to  be  read  into  a  working  depth  matrix  (ZB) . 

The  values  for  KMX  and  KMY,  recorded  cell  resolution 

(DHN) ,  and  the  limits  of  the  input  matrix  ultimately 
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determine  the  axis  ranges.   The  user  should  be  familiar 
with  these  parameters  when  choosing  the  origin  of  the 
plot. 

ZST   -  initial  depth  value  for  vertical  plot  (kilometers) 
ZEND  -  final  depth  value  for  vertical  plot  (kilometers) 
ARST  -  initial  range  value  for  vertical  plot  (kilometers) 
AREND  -  final  range  value  for  vertical  plot  (kilometers) 
The  above  values  are  used  to  draw  the  axes  of  the  vertical 
plot.   The  initial  values  may  be  negative.   The  program 
adjusts  the  final  values  if  necessary  for  plot  esthetics. 
Card  4  (6F8.2) 

HDST   -  heading  of  the  first  ray  in  each  ray  fan  (degrees) 
HDEND  -  heading  of  the  last  ray  in  each  ray  fan  (degrees) 
ELST   -  elevation  angle  of  the  first  ray  fan  (degrees) 
ELEND  -  elevation  angle  of  the  last  ray  fan  (degrees) 
Headings  are  measured  from  grid  north.   Elevations  are 
measured  from  the  horizontal  plane,  positive  downward 
(toward  increasing  depth) . 

HORLIM  -  maximum  ducting  or  channeling  range  (kilometers) . 
This  is  a  gate  value,  like  BLMAX  CARD  2,  which 
terminates  the  trace  of  a  ray  when  exceeded. 
HORLIM  serves  to  screen  out  those  rays  which 
do  not  appear  to  contact  the  bottom. 
SMNGL  -  (optional)  small  angle  value  for  alternative 
cosine  difference  calculation  (degrees) .   A 
value  of  less  than  three  degrees  is  recommended. 
See  the  discussion  on  "small  angle  approximations" 
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in  Chapter  III.   The  default  value  of  two 

degrees  is  set  when  left  blank. 
Card  5  (615) 

NHD    -  number  of  rays  per  ray  fan  (kl) 

NEL    -  number  of  ray  fans  per  bundle (>1) 

NSTOP  -  maximum  number  of  bottom  bounces  per  ray  (>1) 

The  program  restricts  NHD  to  30  rays,  with  a  maximum 

total  number  of  100  rays  produced,  each  with  a  maximum 

number  of  50  bottom  bounces  or  200  reversals. 

INCPLT  -  number  of  ray  segments  defining  one  plot  seg- 
ment (>1) .   With  sufficiently  small  ray  seg- 
ments it  is  not  necessary  to  use  each  iteration 
for  plot  definition  to  achieve  a  smooth  plot, 
so  this  parameter  is  provided  to  enable  the 
user  to  coordinate  FDL  and  INCPLT  for  run 
time  economy. 
NLEN    -  length  of  the  longest  horizontal  plot  axis 

(<19  inches).   This  value  controls  the  physi- 
cal size  of  the  plot.   A  value  less  than  four 
will  supress  the  plot  (resulting  in  an  ABEND 
code  in  the  HASP  log) ,  but  it  will  not  affect 
the  calculations  or  printed  and  punched  out- 
put.  A  value  greater  than  15  requires  addi- 
tional JCL  for  strip  plotting,  since  the 
total  plot  width  is  NLEN  plus  six  inches. 
Refer  to  Ref.  7  for  strip  plotting. 
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NSSP  -  number  of  sound  speed  profiles  (water  masses 
provided  as  data  (  1  <  NSSP  <  3) . 

•  Card  6  (A4,I3) 

OMAS  (1)  -  name  of  first  water  mass 

NVAL  (1)  -  number  of  (depth, sound  speed)  pairs  defining 
the  sound  speed  profile  (1  £  NVAL  £  50) . 

•  Card  7  -  (through  card  16  if  necessary)  (10F8.2) 

These  are  the  data  pairs  (depth,  sound  speed) ,  in  meters 
an  m/sec,  from  the  sea  surface  to  the  last  data  point 
in  the  profile.   The  gradients  will  be  calculated  linearly 
between  these  points,  with  the  final  gradient  below  the 
last  point  set  to  0.01668  m/sec.  by  the  program. 

•  Next  card  (4F8.2) 

If  the  previous  data  is  to  be  followed  by  the  data  for 
another  water  mass,  this  card  must  be  included  to  define 
the  frontal  boundary  between  the  two.   The  order  in 
which  the  water  masses  appear  as  inptut  determines  their 
positions  relative  to  the  fronts.   Cards  which  describe 
the  next  water  mass,  follow  the  frontal  definition  card 
and  are  formatted  as  described  cards  6  and  beyond.   The 
last  data  card  is  the  last  string  of  (depth,  sound  speed) 
pairs  of  the  last  water  mass;  refer  to  Appendix  M  (sec- 
ond page)  for  an  example  of  three  water  masses.   The 
following  variables  define  the  boundary  between  each 
water  mass : 

FENTX  (1)  -  x  coordinate  of  front  beginning  point  (km) 
FENTY  (1)  -  y  coordinate  of  front  beginning  point  (km) 
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FENTX  (1)  -  x  coordinate  of  front  end  point  (km) 
FEXTY  (1)  -  y  coordinate  of  front  end  point  (km) 
The  examples  shown  in  Fig.  IV. 1  identify  the  water 
mass  boundaries  for  various  input  axis  values  in  the 
form  (FENTX,  FENTY) ,  (FENTX,  FEXTY)  for  each  water  mass 
boundary.   The  first  example  a.  represents  the  data 
from  an  original  ECTRACE  run,  refer  to  the  second  page 
of  Appendix  M.   Exambles  b.  and  d.  show  the  priority 
of  the  water  masses.   The  first  water  mass  defined, 
Northwestern  Mohns  Ridge  (NWMR) ,  dominates  the  second, 
Southeastern  Mohns  Ridge  (SEMR) ,  and  the  third,  isove- 
locity  (ISOV) ,  where  they  overlap. 
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FIGURE  IV. 1. 


Examples  of  water  mass  boundaries  for  various  values  of 
(FENTX,FENTY)  and  (FEXTX,FEXTY) . 
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V.   INTERPRETATION  OF  ECTRACE  OUTPUT 

The  following  discussion  refers  to  an  ECTRACE  run  titled 
"Trial  67***Synthetic  Basin***"  (Figure  5) ,  which  is  a  trace 
conducted  in  a  synthetic  parabolic  basin  (Figure  6) .   A  por- 
tion of  the  printed  output  is  shown  in  Tables  II  -  V.   In  this 
ECTRACE  run,  eight  ray  fans  were  selected  to  be  projected  be- 
tween -30  and  +30  degrees  elevation.   Each  ray  fan  was  com- 
posed of  only  one  ray,  thus  the  initial  heading  for  each  ray 
was  010  degrees. 

A.   PRINTED  OUTPUT 

Refer  to  Table  II  through  V  (Appendix  G)  as  examples  of 
the  following  discussion. 

1 .  Parameter  Tableau    (Table  II) 

After  processing  and  adjusting  the  input  parameters 
just  described,  ECTRACE  prints  a  final  tableau  of  these  parame- 
ters plus  some  additional  parameters  calculated  from  the  origi- 
nal input  such  as  layer  gradients,  the  total  number  of  rays 
to  be  traced,  and  the  spacing  between  rays.   In  addition,  the 
results  of  the  depth  matrix-grid  load  is  summarized,  revealing 
data  from  the  header  record  of  the  grid  file. 

2.  Ray  History  Tableau   (Table  III) 

During  the  trace  of  each  individual  ray,  this  tableau 
prints  a  line  of  ray  parameters  at  each  end  point  and  reversal, 
The  following  is  a  description  of  one  line  of  parameters. 
*  TYPE  identifies  the  nature  of  the  reversal,  and  is 
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abbreviated  INIT  for  initial  position,  BOTT  for  bottom 
bounce,  SURF  for  surface  reflection,  OVER  for  turnover, 
and  STOP  for  terminated  position. 

•  X,  Y,  and  DEPTH  are  the  coordinates  in  kilometers  of  the 
reversal  and  end  point. 

•  GRID  HEADING  is  the  ray  heading  measured  clockwise  from 
grid  North  at  the  completion  of  the  reversal. 

•  ELEV  is  the  ray's  elevation  angle  measured  downward 
(positive  toward  increasing  depth)  from  the  horizontal 
zero  degree  reference) . 

•  BRNG  TO  REC  is  the  bearing  from  the  reversal  point  to 
the  receiver  position. 

•  DIST  TO  REC  is  the  slant  range  distance  from  the  reversal 
point  to  the  receiver. 

•  WTR  MASS  is  the  water  mass  identifier  indicating  the 
sound  speed  profile  applied  in  the  present  calculations. 

•  LYR  is  the  constant  gradient  layer  number  corresponding 
to  the  order  input  as  profile  data  (i.e.,  layer  1  is  the 
surface  layer  and  indicates  that  the  surface  gradient 
applies . 

•  DEPTH  LIMITS  are  the  physical  and  refractive  sound  chan- 
nel limits  of  the  ray,  determined  by  its  invariant.   If 
MIN  is  0.0,  an  upward  traveling  ray  will  reflect  off 
the  surface  before  reaching  its  upper  vertex  (turnover) 
depth.   If  the  MAX  depth  is  greater  than  or  equal  to  10km, 
a  line  of  asterisks  is  printed  since  the  ray  could  never 
attain  its  turnunder  depth  before  reaching  the  bottom. 


53 


Otherwise,  these  limits  are  the  ray's  vertex  depths, 
and  if  the  reversal  type  is  OVER  or  UNDR,  the  appropri- 
ate limit  should  exactly  match  the  depth  value  on  the 
same  line. 

•  C  is  the  sound  speed  at  the  reversal  point.   This  is  the 
ray  invariant  used  for  tracing. 

•  GRAZ  ANGLE  is  the  angle  in  degrees  between  the  ray  vector 
and  the  line  of  intersection  of  the  bottom  facet  plane 
with  the  vertical  plane.   This  is  the  angle  used  to 
calcualte  bottom  loss  and  the  new  ray  direction  vector. 

•  BTM  ANGLE  is  the  angle  between  the  bottom  facet  plane 
normal  and  the  vertical  (z-axis) ,  in  degrees. 

•  BTM  LOSS  is  the  accumulated  intensity  loss  in  dB  after 
the  bottom  reflection.   When  this  value  equals  or  exceeds 
BLMAX,  the  ray's  trace  will  terminate. 

•  LPS  is  the  number  of  steps  required  by  the  CONTAC  sub- 
routine to  locate  the  point  of  bottom  contact.   A  long 
ray  segment  length  may  require  more  steps  than  a  shorter 
one.   A  maximum  number  of  50  steps  is  permitted  before 
the  bottom  facet  is  assumed  anomalous  and  the  ray's 
trace  is  terminated. 

•  Prior  to  the  final  line  of  an  individual  ray's  trace 
history,  is  a  summary  of  plot  definition  parameters  plus 
data  on  the  ray's  CPA  to  the  receiver,  and  reason  for 
ray  trace  termination.   The  CPA  calcualtions  are  accurate 
to  within  one  ray  segment  length. 
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3 .  Ray  Fan  Summary  Tableau  (Table  IV) 

At  the  completion  of  the  trace  of  one  complete  ray 
fan,  a  table  is  printed  listing  each  ray's  elevation  angle 
and  depth  in  meters  before  each  bottom  bounce.   All  zero  values 
at  the  bottom  of  each  ray's  column  indicate  blank  data. 

4.  CPA  Summary  Tableau  (Table  V) 

When  all  ray  traces  have  been  completed,  a  table  is 
printed  listing  each  ray's  heading,  time,  and  distance  at  CPA 
to  the  receiver. 

B.   PLOTTED  OUTPUT 

ECTRACE  produces  a  single  montage  of  up  to  five  separate 

plots  .CF.ig..  "13  ,  the  total  dimensions  of  which  depend  upon  the 

user's  selection  for  the  maximum  size  in  inches  of  the  longest 

axis  of  the  horizontal  plot  (NLEN) .   The  following  table  lists 

some  sample  values  for  NLEN  and  the  resulting  dimensions  of 

the  complete  plot. 

Maximum  Maximum 
NLEN   (in.) Height    Width  Remarks 

5  7        8.33      minimum 

10  14       16.67 

15  21       25        maximum  size  without  strip 

plotting  [5] 

19  26.6     31.67     maximum  allowable 

Components  of  the  ECTRACE  plot  follow. 
1.   Sound  Speed  Profiles  (Figure  13) 

A  plot  of  sound  speed  versus  depth  is  made  for  each 
water  mass  defined. 
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2.  Vertical  Plot  (Figs.  15  and  16) 

Most  two-dimensional  ray  trace  programs  produce  plots 
of  ray  paths  on  graphs  of  depth  versus  range.   The  vertical 
plot  is  similar  except  that  all  ray  paths  are  projected  from 
three  dimensions  onto  the  vertical  plane.   The  plane  is  oriented 
along  the  mean  ray  fan  heading  and  includes  its  profile  of  the 
sea  bed.   In  ECTRACE  trial  65B  (Figure  13)  the  ray  fans  each 
contain  two  rays  whose  initial  headings  vary  from  360  to  020 
degrees  as  indicated  by  10.000  +/-  10.000  on  the  plot.   The 
bathymetry  profile  is  made  along  the  heading  of  010  degrees. 
On  a  plot,  the  heading  value  (HEADINGS)  always  indicates  the 
direction  of  the  bathymetry  slice.   The  (+/-)  value  signifies 
the  initial  heading  of  the  first  and  last  rays  and  that  addi- 
tional rays  are  evenly  spaced  in  the  interval. 

Since  the  rays  are  actually  traced  in  three  dimensions, 
bottom  reflections  occurring  outside  the  vertical  plot  plane 
may  appear  to  be  taking  place  above  or  below  the  bottom  sur- 
face instead  of  at  the  boundary  (Fig.  15).   Only  those  bottom 
contact  points  which  occur  in  the  projection  plane  must  lie 
on  the  seabed  profile. 

The  symbols  listed  on  the  plot  (Fig.  5)  are  used  to 
identify  bottom  contact  locations  in  the  vertical  and  horizon- 
tal plot,  and  are  unique  for  each  ray  fan  except  that  up  to 
thirteen  symbols  are  used  before  repetition,  as  shown  in 
Table  VI. 

3.  Horizontal  Plot  (Fig.  17) 

This  plot  is  the  most  important  plot  product  of 
ECTRACE  and  is  the  plot  of  the  ray  paths  on  a  horizontal  plane 

56 


(topview) .   This  plot  displays  the  X-Y  positions  of  each 
bottom  contact  point   using  the  ray's  (ray  fan)  identification 
symbol  and  connects  them  with  straight  lines,  revealing  the 
horizontal  deflection  effect.   Each  ray's  plot  begins  at  the 
source,  ends  up  at  the  last  bottom  contact  point  in  the  grid 
and  is  labeled  with  the  ray  number  (Fig.  17) .   The  ray  num- 
ber can  be  used  to  locate  its  printed  and  punched  output 
history. 

The  axes  are  labeled  with  respect  to  the  origin  of 
the  input  matrix.   The  geographic  center  position  and  coor- 
dinates of  the  input  matrix,  if  included  in  the  bounds  of 
the  working  matrix  (FORTRAN  array  ZB) ,  is  also  plotted  for 
reference.   It  is  important  to  note  that  this  plot,  without 
the  bottom  contact  symbols,  may  combine  the  results  of  other 
ECTRACE  runs  by  combining  the  optional  punched  card  outputs 
as  input  to  ECCOM  (Figs.  7  and  18). 

C.   PUNCHED  OUTPUT 

The  punched  card  output  for  a  single  ECTRACE  job  consists 

of  the  following  cards: 

CARD  DATA  CONTENTS 

1  The  runtitle 

2-3  The  grid  title 

4  Grid  reference  coordinates  and  scale 

factor 

5  Initial  elevation,  initial  heading, 

number  of  points  (BB)  and  ray  num- 
ber of  the  first  ray 

6+  X-Y  coordinate  pairs  of  plot  points 

of  the  first  ray 
(Data  for  additional  rays  is  repeated  as  5  and  6+  above) 
Last  9999.,  to  flag  the  end  of  data 
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These  cards  may  be  stacked  as  a  member  with  those  of 
several  ECTRACE  jobs  as  input  to  a  single  ECCOM  job.   For 
ECCOM  to  be  meaningful,  cards  2  through  4  of  each  stack  mem- 
ber should  be  identical,  and  the  user  should  verify  this 
himself. 

Punched  card  output  may  be  quite  voluminous  and  should 
not  be  produced  unless  required.   The  following  card,  inserted 
before  the  first  GO  card  (JCL) ,  will  suppress  the  punched 
output: 

//GO.FT07F0  01  DD  DUMMY 
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VI.       SUMMARY 

A.  POTENTIAL  APPLICATIONS 

ECTRACE,  in  its  present  form,  has  many  potential 
applications : 

•  The  investigation  of  asympotic  deflection  angles  in 
ocean  areas  where  bathymetry  data  are  well-documented, 
allows  optimum  hydrophone  positioning.   For  example, 
Ref.  6  cites  cases  where  seemingly  attractive  receiver 
positions  would  be  affected  by  shadowing  in  the  hori- 
zontal plane. 

•  Adaptation  of  ECTRACE  ray  tracing  methods  could  be 
incorporated  into  existing  programs  which  presently 
rely  upon  an  assumption  of  radial  symmetry  of  ray 
propagation  about  a  source  or  receiver  (no  horizontal 
deflection) . 

•  The  horizontal  deflection  effect  may  lend  additional 
insight  into  the  study  of  travel  times  and  intensities 
of  acoustic  signals,  enabling  more  accurate  source 
fixing  information. 

B.  IMPROVEMENT 

There  are  many  areas  of  improvement  which  would  increase 
the  utility  of  ECTRACE  and  its  supporting  computer  programs 
for  the  above  applications: 

•  A  better  routine  to  convert  bathymetry  contour  data 
into  a  matrix  of  depth  values  is  essential  for  modeling 
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regions  for  which  only  contour  data  exists.   Most  known 
interpolation  algorithms  perform  best  on  data  which  are 
randomly  or  regularly  spaced  rather  than  in  the  form  of 
contours.   For  example  RAIN 2  [5]  would  probably  perform 
better  with  the  data  used  to  make  the  contours  than  it 
does  with  the  contours  themselves.   While  the  contour 
data  are  useful  for  cartographic  and  geologic  purposes, 
bathymetry  data  in  a  less  processed  form  are  significantly 
better  sources  for  constructing  computer  models  of  the 
sea  bed. 

ECTRACE  in  its  present  form  makes  no  calculations  of 
intensities  or  their  losses  due  to  spreading,  attenua- 
tion, or  scattering,  thus  preventing  its  use  as  a 
propagation  loss  model. 

An  improved  bottom  loss  function  is  recommended  over 
the  supplied  FORTRAN  function  subprogram  FBLOSS.   The 
user  may  perform  this  substitution  easily  by  deleting 
FBLOSS  from  the  INCLUDE  step  in  the  LINK  statement  of 
the  ECTRACE  job  control  and  supplying  another  function 
subprogram  (added  at  the  end  of  the  calling  program) 
that  accepts  a  double-precision  argument  for  grazing 
angle  and  is  also  called  FBLOSS. 

The  assumption  of  specular  reflection  neglects  the 
phenomenon  of  transmission  into  the  sediment  and  its 
effect  on  the  location  of  the  reflected  ray  due  to  sub- 
bottom  boundaries.   It  would  be  interesting  to  incorporate 
in  an  improved  version  of  ECTRACE. 
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APPENDIX  A 
GENBOT  USER  INSTRUCTIONS  AND  OUTPUT  DESCRIPTION 

GENBOT  is  used  in  the  construction  of  a  depth  grid  file  from  bathymetry 
contour  files.   In  this  discussion,  the  term  "grid"  refers  to  a  matrix  of 
depth  values  whose  row  and  column  dimensions  represent  data  point  spacing 
in  terms  of  a  desired  resolution  parameter  DHN. 
1.   Selection  of  Bathymetry  Grid  Regions 

GENBOT  will  create  a  bathymetry  grid  file  for  any  region  of  the  geoid 
for  which  contour  data  is  supplied  in  the  format  stated  in  section  A-4. 
The  user  must  select  the  latitude  and  longitude  of  the  center  of  the 
region  of  interest,  a  radius  in  kilometers  and  the  cell  resolution  in 
kilometers.   GENBOT  will  attempt  to  create  a  square  matrix  whose  semi- 
axes  represent  the  selected  radius  and  whose  elements  represent  depth 
values  in  kilometers. 

A  30-km  data  margin,  outside  the  calculated  axis  limits,  is  required 
to  begin  interpolation.   Should  this  requirement  not  be  met,  GENBOT  will 
symmetrically  reduce  the  length  of  one  or  both  axes  until  the  resulting 
rectangular  matrix  plus  margins  with  the  desired  center  position  intact 
is  totally  encompassed  within  the  data  region.   The  center  coordinates, 
area  radius,  and  cell  resolution  must  be  selected  so  as  to  avoid  generating 
a  grid  which  is  impracticably  small.   In  particular,  a  generated  grid  may 
not  overlap  a  pole,  and  one  with  a  pole  on  a  margin  boundary  must  be  sup- 
plied data  from  eighteen  contour  files  to  cover  180  degrees  of  longitude. 

Contour  data  currently  available  in  usable  format  cover  40  degrees 
of  longitude  of  the  North  Atlantic  Ocean,  specifically  the  Norwegian  and 
Greenland  Seas,  from  longitudes  20W  to  20E  and  latitudes  60N  to  90N 
(Fig.  8).   The  Naval  Research  Laboratory,  Acoustics  Division,  Arctic  Section, 
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has  additional  files  which  cover  the  entire  North  Atlantic  land  and  sea 
regions  which  may  be  utilized  by  GENBOT  after  modification.   However, 
GENBOT  need  not  be  the  exclusive  means  of  grid  instruction.   The  only 
requirement  is  that  the  grid  exists  as  a  file  in  the  format  specified. 
The  user  may  wish  to  draw  on  other  souces  of  bathymetry  data,  such  as 
SYNBAPS  [8],  and  tailor  other  available  resources  to  construct  a  usable 
grid. 
2.   Grid  Construction  Techniques 

Once  supplied  with  the  user's  selection  for  the  grid  title,  center 
coordinates,  cell  resolution  and  area  radius  (TTL,  CNLAT,  CNLON,  DHN, 
ARAD) ,  GENBOT  begins  by  defining  a  north-oriented  square  region  whose 
boundaries  are  at  a  distance  ARAD  from  the  center.   The  subscripts  of 
the  two  dimension  matrix  depicting  this  region  will  represent  equal  units 
of  DHN  spacing  from  the  origin  along  the  x  and  y  axes  of  the  square  using 
a  cartesian  coordinate  system.   Thus,  while  all  distances  calculated 
during  grid  generations  are  geodetic  distances  from  the  center,  point 
positions  are  made  relative  to  the  grid  origin  or  southwest  corner.   This 
technique  results  in  a  plane  projection  which  is  distortion  free  along 
straight  tracks  which  pass  through  the  center.   Subroutine  CORNER  then 
calculates  and  prints  a  table  of  the  initial  cartesian  and  geographic 
coordinates  of  the  grid  corners.   The  cartesian  coordinates  will  always 
range  from  (0,0)  to  their  maxium  values,  initially  double  the  value  of 
ARAD.   The  geographic  coordinates  should  reflect  some  meridian  conver- 
gency.   This  table  will  be  repeated  for  the  final  grid,  showing  the 
adjustments  made,  if  any,  to  meet  the  data  requirements. 

Many  geographic  coordinate  pairs  in  an  input  file  may  be  well  outside 
the  area  relevant  to  the  grid,  so  GENBOT  calculates  gate  values  called 
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SLAT,  ELAT,  SLON,  and  ELON  for  data  discrimination.   SLAT  and  ELAT  are 
always  one  degree  beyond  the  minimum  and  maximum  latitude  values  calcu- 
lated by  CORNER,  but  SLON  and  ELON  may  extend  beyond  the  minimum  and 
maximum  longitudes  by  1,  10,  or  20  degrees,  depending  upon  the  maximum 
latitude. 

The  grid  region  and  its  data  environs  recognized  by  GENBOT  at  this 
stage  is  shown  in  Fig.  A.l.  All  points  within  these  environs  are  sub- 
ject for  consideration  during  interpolation  of  depth  values  in  the  grid 


FIGURE  A.l.   Geographic  region  used  for  grid  construction, 


63 


Data  files  used  as  input  for  GENBOT  are  digitized  represenations  of 
an  NRL  bathymetry  chart  [4]  and  were  recorded  on  Mobil  Oil  Company  Tape 
NAV-78  for  use  in  a  computer  program  to  draw  portions  of  the  chart  in  the 
form  of  bathymetric  contour  lines  and  their  depth  labels  as  they  exist 
on  the  chart.   The  records  are  in  the  card  image  format  and  consist  of 
three  types: 

•  NAME  -  depth  value  coordinates 

•  CON  -  contour  line  header  label 

•  Chart  coordinates  of  points  along  the  contour  line 

In  the  original  form  the  numerical  data  on  the  tape  consisted  of  lati- 
tude/longitude pairs,  providing  no  means  of  determining  the  depth  repre- 
sented by  the  contour  lines  without  visual  inspection  of  the  computer  plot 
The  tape  used  in  this  research  is  similar  in  format  to  NAV-78  except  that 
the  contour  line  header  records  include  depth  data  (provided  by  NRL) ,  and 
only  the  regions  shown  in  Fig.  A. 2  and  listed  in  Table  A.I,  are  covered. 
The  numerals  of  the  data  set  names  were  chosen  to  reflect  the  original 
file  number  on  the  source  tape,  NAV-78. 

Using  the  magnetic  tape  files  selected  by  the  user,  GENBOT  begins 
reading  the  digitized  contour  data,  stopping  after  each  record  to  dis- 
criminate the  points,  calculate  geodetic  distances,  transform  the  data 
to  x,  y,  z  coordinates,  and  load  the  coordinates  into  three  axis  vectors, 
AX,  AY,  and  AZ  (which  will  be  used  later  for  interpolation)  while  simul- 
taneously keeping  a  running  plot  of  the  contour  lines. 

After  all  input  files  have  been  processed,  the  vectors  AX  and  AY  are 
sorted  into  ascending  order  by  the  SSP3  subroutine  RAIN1  [5].  The  mini- 
mum and  maximum  values  are  tested  to  determine  if  the  30-km  outside  data 
margin  is  included  in  the  range.   If  not,  axis  reduction  begins  until 
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72' 


70° 


68° 


FIGURE  A. 2.   North  Atlantic  regions,  with  their  magnetic  tape  file  numbers 

and  data  set  names,  available  for  use  by  GENBOT.   All  are  10° 

by  10° .   Other  regions  available  on  the  same  tape  are  listed 
in  Table  A.I. 
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TABLE  A.I 
FILE  ONE  OF  CONTOUR  DATA  RECORDED  ON  A  MAGNETIC  TAPE 


DATA:     BATHYMETRIC  CONTOURS  OF  THE  NORWEGIAN  AND  GREENLAND  SEAS  FROM 
LAT  60N  TO  90N  and  LONG  20W  to  20E. 

SOURCE:   ADAPTED  FROM  MOBILE  OIL  CO.  TAPE  NO.  NAV-78  AND  NAVAL  RESEARCH 
LABORATORY  DATA. 

STRUCTURE:   EACH  FILE  COVERS  A  10  X  10  DEGREE  REGION.   THE  FOLLOWING  LIST 

IDENTIFIES  THESE  FILES.   THE  FIRST  VALUE  IS  THE  FILE  NUMBER,  FOLLOWED  BY 

THE  DSN,  THEN  THE  LAT /LONG  PAIR  OF  THE  SOUTHWEST  CORNER  (LABOR, LONOR)  OF 

THE  REGION,  AND  FINALLY  THE  NUMBER  OF  CARD  IMAGES  CONTAINED  IN  THE  FILE. 

THE  TAPE  IS  STANDARD  LABELED,  DENSITY  1600  BPI ,  LRECL=80 ,BLKSIZE=800. 


2 

CONF13 

+80-20 

256 

5 

CONF26 

+80-20 

4343 

8 

CONF34 

+70+00 

12443 

11 

CONF35 

+70-10 

9428 

3 

CONF24 

+60+10 

569 

6 

CONF27 

+60-10 

6001 

9 

CONF37 

+80+00 

491 

12 

CONF36 

+80+10 

409 

4 

CONF25 

+60+00 

3026 

7 

CONF30 

+70+10 

857 

10 

CONF06 

+70-20 

4454 

13 

CONF38 

+80-10 

972 
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either  the  margin  requirement  is  met  or  goes  below  a  7.5  km  semi-axis 
length,  the  latter  case  preventing  grid  construction  although  a  contour 
plot  may  still  be  produced.   Subroutine  CORNER  is  then  called  again  to 
print  the  corner  coordinate  table.   If  no  adjustments  were  necessary,  the 
second  table  will  be  identical  to  the  first. 

The  interpolation  process  begins  at  the  southwest  corner  of  the  grid 
and  proceeds  in  a  columnar  direction  (i.e.,  completing  the  range  of  x. 
values  for  each  y.  value)  via  successive  calls  to  the  SSP3  subroutine 
RAIN2.   The  accuracy  of  the  interpolated  values  is  proportional  to  the 
density  of  data  in  the  given  area  of  interest.   Every  ten  values  inter- 
polated are  written  as  a  record  on  the  output  file.   The  process  continues 
until  the  total  number  of  matrix  elements  has  been  reached,  at  which  time 
the  output  file  is  closed  and  the  program  is  terminated. 
3.   User  Instructions 

a.   JCL  and  Program  Deck. 

The  program  section  of  a  sample  deck  is  printed  below.   The  TIME 
and  REGION. GO  parameters  were  chosen  in  anticipation  of  a  maximum  number 
of  5000  data  points  in  the  area  of  interest,  a  square  region  of  210  km 
on  a  side  (150  km  plus  margin) .   The  program  listings  included  in  this 
appendix  are  part  of  this  example. 

//(Standard  job  card) ,TIME=20 

//  EXEC  FORTCLGW, REGION. GO=250K  (See  Appendix  0.) 

//  FORT. SYS IN  DD  * 

(CENBOT,  CORNER,  GEOPLT,  GEODST  listings) 
(RAIN1,  RAIN3,  RAIN2  listings) 

/* 

The  following  DD  statements  identify  the  magnetic  tape  files 
representing  the  ten  by  ten  degree  regions  into  which  the  subroutine 
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extends.   In  the  example,  the  area  chosen  is  centered  at  coordinates 
70°-04,N/15o-51'  (70.07059/15.85302)  with  a  radius  of  75  km.   This  area 
requires  two  data  files  as  input,  CONF24  (file  3)  and  CONF30  (file  7). 
Since  tape  data  is  sequential,  the  order  of  the  DD  cards  should  be  in  the 
order  of  the  file  numbers  (appearing  in  the  LABEL  parameter) . 
//GO.FT01F001  DD  UNIT=3400-3,  V0L=SER=NPS,DISP=(0LD,PASS> , 
//   LABEL=(3,SL,,IN) ,DSN=CONF30 

//GO.FT01F002  UNIT=AFF=FT01F001,VOL=(, RETAIN, REF=*. GO. FT01F001) , 
//    DISP=(0LD,PASS) ,LABEL=(7,SL, ,IN) ,DSN=C0NF30 
More  information  on  magnetic  tape  usage  at  NPS  may  be  found  in  Ref .  10. 
Next  must  follow  the  DD  statements  specifying  the  data  set  and 
device  for  permanent  storage  of  the  generated  grid.   The  permanent  device 
may  be  another  magnetic  tape,  a  disk,  or  punched  cards.   In  all  cases, 
the  size  of  the  output  data  set  is  important,  and  the  DCB  parameters 
should  be  considered  in  advance  if  the  grid  size  is  different  from  the 
example.   GENBOT  produces  card  image  out  put,  meaning  the  logical  record 
length  is  80  bytes. 

•  Punched  cards.   If  the  grid  is  to  be  stored  on  punched  cards,  the 
user  need  only  change  the  data  set  reference  number  on  the  appro- 
priate WRITE  statements  in  GENBOT  from  a  two  to  a  seven.   The  num- 
ber of  cards  punched  will  equal  a  3+[(row  dimension) (column 
dimension) /10]  rounded  up  to  an  integer  value.   This  represents 
three  header  cards  followed  by  ten  depth  values  per  card  in  columnar 
order,  as  required  by  ECTRACE.   Thus  if  the  input  grid  radius  was 
75  km  and  input  cell  resolution  1  km,  the  maximum  parameter  calculated 

by  GENBOT  will  be: 

XMAX  YMAX  IMAX  JMAX 

150  150  151  151 
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resulting  in  a  matrix  of  22801  elements  and  a  punched  output  of 
2284  cards. 

•  Disk.   Grid  storage  on  a  disk  device  is  the  most  convenient  form  for 
multiple  ECTRACE  runs  over  the  same  sea  bed.   The  following  DD  state- 
ment stores  the  example  grid  on  a  3330  disk  drive  device  at  NPS  as  a 
data  set  named  S999.GRID. 

//GO.FT02F001  DD  UNIT=3330,VOL=SER=DISK03,DISP=(,KEEP) , 

//   DCB=(FECFM=FB,LRECL=80,BLKSIZE=800) ,SPACE= (TRK, 20) , 

//   DSN=S9999.GRID 

Blocked  at  800  bytes,  the  output  data  requires  14  tracks  per  10 
records  (card  images) .   More  information  on  the  DCB  and  DSN  parameters 
may  be  found  in  Ref.  11.   When  the  generated  grid  is  no  longer  needed, 
and  if  the  user  wishes  to  generate  a  new  grid,  or  if  space  has  already 
been  allocated  on  the  device,  the  DISP  parameter  should  be  SHR  in 
place  of  (KEEP)  and  the  DCB  and  SPACE  specifications  may  be  omitted. 
Disk  data  sets  are  periodically  purged  by  center  operations  personnel, 
but  obsolete  data  set  space  should  be  cleared  or  reallocated  by  the 
user. 

b.   Input  Data. 

The  following  input  data  are  user  selected  values  and  follow 
the  //GO.SYSIN  DD*  card. 


Card 

Variable 

FORTRAN 

Field  Descriptors 

1-2 

TTL 

6A8/6A8 

3 

CNLAT , 
CNLON , 
DHN,  ARAD 

4F10.5 

NFILES,  LEN 


215 


Remarks 

Title 

center  latitude, 
center  longitude, 
cell  resolution, 
radius 

Number  of  input 
files,  plot 
length 
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The  magnitude  of  LEN  determines  the  plot  size  (.<19  inches)  .  If 
LEN  is  negative,  no  grid  will  be  constructed.  If  its  magnitude  is  less 
than  four,  no  plot  will  be  produced. 

c.  Printed  Output 

GENBOT  printed  output  supplies  the  following: 

•  A  repeat  of  the  user  input  data. 

•  A  table  of  the  grid  corner  coordinates  and  the  depth  matrix 
dimensions  using  the  input  parameters. 

•  A  summary  of  results  for  each  contour  line,  including  its  depth, 
and  the  title  from  the  header  label,  the  number  of  points  on  the 
line  (NPTS)  and  the  number  of  points  used  from  the  line  (KPTS) . 

•  A  corner  coordinate  table  of  the  final  grid. 

•  A  summary  of  the  interpolation  results.  A  line  is  printed  for 
each  interpolated  point  whose  accuracy  may  be  questioned.  The 
meaning  of  the  IER  parameter  may  be  found  in  the  RAINl  listing 
in  Appendix  I. 

d.  Plotted  Output 

The  contour  lines  of  the  input  file  are  plotted  with  the  local 
parallels  and  meridians  overlaid.   The  resultant  subregion  represented  by 
the  output  grid  is  outlined  by  a  rectangle  with  labeled  axes.   The  axis 
values  are  referenced  to  the  desired  grid  origin  and  not  the  resultant 
grid  origin,  the  latter  being  at  (0,0)  always. 

e.  Programming  Notes  on  Contour  Data  Interpolation 

Subroutine  RAIN2,  the  interpolation  routine  used  by  GENBOT,  is 
aptly  suited  to  interpolate  points  that  are  more  or  less  randomly  spread 
in  both  horizontal  directions.   The  nature  of  contour  data,  however,  is 
typically  a  collection  of  closely  spaced  depth  values  along  widely  spaced 
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contour  (compare  Figures  III . 1.  and  A. 1. ) .   This  type  of  data  arrangement 
frequently  results  in  a  colinearity  problem  for  RAIN2,  which  can  be  resolved 
only  by  searching  ever  widening  data  cells  for  a  non-colinear  point,  result- 
ing in  a  degradation  of  accuracy.   In  addition,  the  routine  is  highly  sensi- 
tive to  contour  line  kinks  which  result  in  right-angle  deviations  in  the 
interpolated  contour  lines  (which  may  be  plotted  by  GRDCHK) . 

Correcting  the  deficiencies  of  an  interpolating  routine  is  beyond 
the  scope  of  this  research,  but  the  Naval  Research  Laboratory,  Arctic  Sec- 
tion, is  currently  investigating  methods  for  generating  grid-matrices  from 
contour  data.   Users  are  encouraged  to  modify  the  code  used  in  GENBOT  as 
desired  to  produce  results  which  more  closely  reflect  the  actual  data 
supplied. 
4.   Input  Contour  Data  Format 

As  stated  previously,  some  contour  data  already  exists  in  a  format 
usable  by  GENBOT,  but  if  data  is  supplied  by  another  source,  it  must  be 
arranged  in  a  file  in  the  manner  and  format  specified  below.   All  records 
are  80  byte  card  images. 

The  first  record  of  the  file  contains  data  on  the  file  itself.   GENBOT 
uses  a  415  field  to  represent  the  origin  latitude/longitude  pair,  the  num- 
ber of  contour  lines  in  the  file,  and  the  number  of  card  images  it  may 
expect  to  read. 

Contour  line  data  exists  as  separate  sequences  of  cards  and  may 
appear  in  any  part  of  the  file  as  long  as  each  line  or  line  segment  is 
in  an  unbroken  sequence  with  a  header  card.   The  header  card  must  contain 
the  line's  depth  in  hundreds  of  meters,  right  justified  in  a  14  field  in 
columns  41-44.   The  cards  following  the  header  card  contain  the  latitude/ 
longitude  pairs  in  degrees  and  decimals  of  degrees  of  the  points  along 
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the  contour  line  in  the  form  of  four  pairs  of  coordinates  (in  FORTRAN 
fields  8F10.5)  per  card.   South  and  west  values  are  negative.   A  value 
greater  than  90.0  in  magnitude  in  any  latitude  field  is  interpreted  as  the 
end  of  data  for  the  contour  line,  and  GENBOT  expects  the  next  record  to 
be  either  another  contour  line  header  label,  the  end  of  the  file,  or  a 
depth  value  position  record. 

Depth  value  position  records  are  optional  and  are  used  during  the 
chart-drawing  phase  of  the  program  to  place  alphanumeric  symbols  on  the 
chart  product.   They  may  exist  any  place  after  the  first  record  in  the 
file  except  imbedded  in  a  contour  line  sequence.   Columns  4-8  must  con- 
tain "NAME"  in  an  A4  field,  columns  19-38  contain  the  latitude  and  longi- 
tude of  the  symbols  in  a  2F10.5  field  and  the  symbols  themselves  must  be 
in  columns  46-47.   These  records  are  normally  used  to  draw  depth  values 
in  hundreds  of  meters  on  the  contour  lines. 

If  the  user  wishes  to  use  GENBOT  on  externally  supplied  contour  data, 
it  it  recommended  that  a  dump  be  performed  on  one  of  the  smaller  files 
(such  as  C0NF13)  of  the  GENBOT  data  tape  to  verify  the  required  format 
for  adaptation. 
5.   Output  Grid  Data  Format 

GENBOT  constructs  the  output  file  in  accordance  with  GRDSCT  require- 
ments.  Since  GRDSCT  is  the  subroutine  used  by  TRACER  (of  ECTRACE), 
G3DP  and  GRDCHK,  any  externally  supplied  grid  data  must  all  be  of  the 
same  structure  as  outlined  below  in  order  to  be  used  by  those  programs. 
Again,  all  records  are  card  images. 

a.  The  first  two  records  contain  the  input  grid  title  in  the  1QA8 

fields. 

b.  The  third  record  contains  the  input  latitude  and  lontitude  of  the 
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grid  center,  the  point  spacing,  and  the  number  of  rows  and  columns  of  the 
matrix  representing  the  grid.   The  record  format  is  (3F10. 5,215) . 

c.   The  fourth  and  remaining  records  contain  ten  depth  values  per 
card,  in  kilometers,  as  10F8.4  fields,  beginning  with  column  one  (y  =  0) 
and  proceeding  all  x-values  before  going  to  the  next  column,  in  accordance 
with  the  FORTRAN  convention  of  columnar  storage  of  matrices. 
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APPENDIX  B 
SYNGEN  DESCRIPTION  AND  USER  INSTRUCTIONS 


1.   Program  Description 

a.  SYNGEN  is  a  collection  of  simple  functions  which  may  be  selectively 
used  to  generate  and  record  a  matrix  of  depth  values  representing  a  desired 
synthetic  bottom  configuration.   Six  bottom  cases  are  provided  (see  below), 
but  the  program  may  be  easily  modified  by  the  user  to  produce  a  case  not 
described  here.   A  program  listing  is  provided  in  Appendix  H. 

The  user  selected  parameters  are: 

•  CNLAT  -  any  desired  center  latitude  value. 

•  CNLON  -  any  desired  center  longitude  value.   These  parameters 
are  requested  in  compliance  wittr  the  required  output  file 
format . 

•  DMIN  -  the  minimum  depth  value  (D  )  to  be  attained  within  a 
75  km  radius  from  the  grid  center. 

•  DMAX  -  the  maximum  depth  value  (Dj  within  the  same  radius. 

•  ITYPE  -  an  integer  code  for  the  bottom  type  desired. 

b.  Table  B.l  lists  the  ITYPE  code,  the  bottom  type,  the  functions 
used  to  generate  the  depth  matrix,  and  the  matrix  dimensions.   All  grids 
represent  a  region  covering  a  150  km  in  both  the  x  and  y  directions,  but 
the  matrix  size  will  be  11  by  11  or  151  by  151,  depending  upon  the  bottom 
type  (DHN  values  are  15  km  and  1  km,  respectively).   In  all  cases, 

AD  =  D  -  D  . 

i     o 

c.  The  depth  file  provided  by  SYNGEN  conforms  to  the  format  required 
by  ECTRACE,  GRDCHK,  and  G3DP.   When  one  of  these  programs  uses  a  SYNGEN 
product,  the  array  dimension  in  the  program's  source  code  must  be  checked 
to  conform  to  the  dimensions  of  the  generated  matrix. 
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TABLE  B.I 
SYNTHETIC  SEA  BED  FUNCTIONS 

Matrix 
ITYPE  Type        Function  Parameters  Dimensions 


12  Wedge      H  =  D  -  y  y     =  AD/135  11  x  11 

0  £  x  £  135 

13  Ridge       H  =  D  +|y$(x)|       y   =  AD/15  11  X  11 

o 

H  =  D  elsewhere     6(x)  =  x  -  75 
o 

60  £  x  <,   90 


14    Trough      H  =  D  -|y5(x)|       6  =  AD/15  11  X  11 

H  =  Di  elsewhwere     <5(x)  =  x  -  75 

60  <  x  <  90 


15    Conical     H  =  D  +  Yr  Y  =  AD/75  151  X  151 

o 


Seamount 


H  =  D  elsewhere     r  =  \Ix  -  75)  +  (y  -75) 


o 


=|7T 


2L     /     *     -,cn2 


16  Sinusoidal     H  =  D  yU+cos{5(x)/R}]  Y  -  Ad/2  151   X  151 

Undulation 

6(x)    =  2(x  -    75) 

R  =  15 


17    Parabolic   H  -■  D_    yr2  Y  =  AD/(75)2  151  X  151 

r  =  >|(x  -  75)  2+  (y  -  75) 2 


basin 
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2.   User  Instructions 

THE  JCL  of  a  sample  deck  are  arranged  as  follows: 

//(Standard  JOB  card) ,TIME=2 

//  EXEC  FORTCLG,  REGION .GO=250K 

//FORT. SYS IN  DD  * 

(SYNGEN  source  code) 

/* 

//GO.FT02F001  DD  UNIT=3330, VOL=SER=DISK03, 

//  DSN=S9999.SYNBED,  DCB= (RECFM=FB,LRECL=80,BLKSIZE=800) , 

//  DSP= (NEW, KEEP),  SPACE= (TRK, 20) 

//GO.SYSIN  DD  * 

(user  selected  data) 

/* 

Notice  that  the  output  data  set  is  directed  to  DISK03.   The  DSN  is  the 
user's  selection  for  the  data  set  name  in  accordance  with  NPS  system  re- 
quirements.  A  151  x  151  depth  matrix  requires  twenty  tracks  of  disk 
space  on  a  3330  device. 

Only  one  data  card  is  needed  to  supply  the  program  with  the  user's 
parameter  selections.   Its  field  descriptor  is  (4F10.5,I5)  and  contains 
the  values  reflected  in  the  same  order  as  listed  in  paragraph  l.a. 
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APPENDIX  C 
G3DP  DESCRIPTION  AND  USER  INSTRUCTIONS 


1.   Program  Description 

a.  G3DP  uses  the  NPS  library  subroutine  PLT3D1  [Ref.  12]  to  draw  a 
perspective  surface  plot  of  a  subset  of  a  depth  matrix.  This  plot  is  a 
normal  projection  of  the  surface  using  a  focal  length  and  line-of-sight 
calculated  from  the  user's  selections  for  the  viewpoint  position  angles 
DEGUP  and  DEGCW.  As  the  variable  names  imply,  the  viewer  is  elevated 
DEGUP  degrees  from  a  reference  depth  and  rotated  DEGCW  clockwise  from  a 
south-to-north  viewing  aspect.  For  example,  the  standard  working  matrix 
(dimensioned  151  by  151)  is  best  viewed  with  DEGUP  and  DEGCW  values  of 

30  degrees  and  45  degrees  respectively,  providing  the  same  viewing  aspect 
demonstrated  in  all  the  the  G3DP  products  printed  in  this  report.   Using 
a  DEGUP  value  of  zero  results  in  an  isometric  rather  than  a  perspective 
projection.   The  DEGCW  selection  should  be  between  zero  and  89  degrees 
(south-north  to  west-east  viewing  aspects) . 

b.  As  in  ECTRACE  and  GRDCHK,  G3DP  requires  source  code  modification 
by  the  user  to  tailor  the  program  to  the  desired  dimensions  of  the  work- 
ing matrix.  Instructions  in  this  regard  are  included  as  comments  in  the 
FORTRAN  listing  (Appendix  K) . 

c.  An  important  characteristic  of  G3DP  is  that  it  uses  only  the 
odd-numbered  rows  and  columns  of  the  working  matrix  to  produce  the  sur- 
face plot.   This  reduction  in  resolution  is  required  to  reduce  core  re- 
quirements and  unclutter  the  graphic  output  caused  by  typically  large 
(greater  than  75  by  75)  working  matrices. 
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2.   User  Instructions 

The  JCL  of  a  sample  program  deck  are  arranged  as  follows: 

//(Standard  JOB  card) ,TIME=2 

//EXEC  FORT CLGW, REGION. G0= 3 00K   (Appendix  0.) 

//FORT. SYS IN  DD* 

(G3DP  source  code) 
/* 

//LINK.USDD  UNIT=3330,VOL=SER=DISK01,DSN=1570.  ECTRACE, 
//   DISP=SHR,LABEL=(,  ,  ,IN) 
//LINK.SYSIN  DD* 

INCLUDE  USDD(GRDSCT) 

ENTRY  MAIN 
/* 
//GO.FT01F 

//GO.SYSIN  DD* 

(user-selected  data) 

/* 

The  LINK  step  provides  access  to  GRDSCT  for  loading  the  working 
matrix  from  the  input  file.   The  input  data  set  is  the  same  example  used 
in  ECTRACE  and  GRDCHK.   The  card  input  data  are  user  selected  and  con- 
tained on  one  card,  as  follows  (FORTRAN  field  descriptor  5F8.3) 

•  XST  -  working  grid  matrix  x-origin  relative  to  that  of  the 
input  grid  (km) . 

•  YST  -  working  grid-matrix  y-origin  relative  to  that  of  the  input 
grid  (km) . 

•  DEGUP  -  viewpoint  elevation  angle  (0  £  DEGUP  <  89  degrees). 
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•  DEGCW  -  viewpoint  rotation  angle  from  south  to  north 
(0  £  DEGCW  <  89  degrees) . 

•  LW  -  width  of  the  plot  (5  <  LW  <  18  inches). 
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APPENDIX  D 
GRDCHK-USER  INSTRUCTIONS  AND  OUTPUT  DESCRIPTION 


1 .  Purpose  and  Output  Description 

Cases  may  arise  where  the  accuracy  of  depth  values  in  a  generated  ma- 
trix may  be  in  question  because  of  data  sparsity,  interpolation  deficiencies 
or  other  reasons.   Program  GRDCHK  assists  in  evaluating  the  generated  grid 
by  printing  and  plotting  selected  columns  of  the  matrix  for  comparison  with 
that  of  the  source  contours.   The  contour  plot  may  also  be  used  on  a  light 
table  with  an  ECTRACE  or  ECCOM  horizontal-ray  trace  plot  (Figs.  5,  6,  and 
7). 

The  NPS  library  routine  CONTUR  is  called  for  producing  the  contour 
plot.   Subroutine  GEOPLR  is  a  modified  version  of  GEOPLT  to  accommodate 
the  rotated  plot  produced  by  CONTUR. . 

2.  User  Instructions 

To  aid  in  this  discussion  refer  to  a  computer  listing  of  GRDCHK  code 
in  Appendix  J.   A  listing  of  the  CONTUR  subroutine  may  be  found  in  Ref.  6. 

The  following  example  applies  to  a  standard  151  by  151  working  grid. 
Below  is  a  listing  of  the  JCL  cards  and  the  position  of  the  program  deck. 

//(Standard  job  card),  TIME=5 

//  EXEC  FORTCGLW, REGION. GO=350K 

//FORT. SYS IN  DD  * 

(GRDCHK, CORNER, GEOPLT, GEODST  source  decks) 

/* 

The  following  DD  statements  link  to  a  member  of  ECTRACE  for  allocating 

a  working  subregion  of  the  input. 

//LINK.USDD  DD  UNIT=3330,VOL=SER=DISK01 , 
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//  DSN=1570.ECTRACE,DISP=SHR,LABEL=(,,,IN) 
//LINK.SYSIN  DD  * 

INCLUDE  USDD  (GRDSCT) 
ENTRY  MAIN 
/* 

Next  follows  the  DD  statement  identifying  the  depth  matrix  input  file, 
The  same  example  is  used  for  the  GENBOT  output  file. 

//GO.FT0F001  DD  DSN=S9999 .GRID,UNIT=3330,VOL=SER=DISK03, 
//    DISP=SHR,  LABEL- (,,, IN) 
Finally  the  user  input  parameters  are  included.   The  card 

//GO. SYS IN  DD  * 
is  followed  by: 

Card  Parameter  Names         FORTRAN  field  description 

1  NL,  LEN  215 

2  ZST,  ZEND,  XMIN,  4F6.1 
YMIN 


where 


3  IPS,  IPE,  JPS,  JPE  415 

•  NL  is  the  number  of  contour  levels  desired.   The  shallowest  and 
deepest  levels,  if  known,  should  be  read  into  ZST  and  ZEND.   If 
the  boundary  levels  are  not  known,  NL  should  be  made  negative, 
and  subroutine  CONTUR  will  calculate  the  boundary  and  contour 
levels.   If  NL  is  greater  in  magnitude  than  100,  no  contour  plot 
will  be  produced. 

•  LEN  is  the  length  of  the  contour  plot's  longest  side  in  inches. 
If  LEN  is  negative,  a  cartesian  grid  will  be  drawn  on  the  con- 
tour plot  instead  of  geographic  coordinate  lines.   If  LEN  is 
less  than  4,  no  contour  plot  will  be  produced.   The  maximum 
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magnitude  of  LEN  is  19  inches. 

•  ZST  AND  ZEND  are  estimates  of  the  smallest  and  largest  values  in 
the  working  matrix  (ZB) . 

•  XMIN  and  YMIN  establish  the  origin  of  the  working  matrix  relative 
to  that  of  the  input  matrix.   The  program  assignments  for  KMX  and 
KMY,  the  recorded  cell  spacing  DHN,  and  the  dimensions  of  the  work- 
ing matrix  determine  the  subregion  size  parameters  in  the  same  man- 
ner as  they  do  in  ECTRACE,  since  subroutine  GRDSCT  is  linked  to 
this  program.   In  effect,  the  user  must  tailor  certain  portions 

of  the  program  cards  to  compensate  for  the  FORTRAN  G  inability 
to  perform  object  time  dimensioning. 

•  IPS,  IPE,  JPS  and  JPE  identify  the  beginning  and  end  rows  and 
columns  of  the  working  depth  matrix  to  be  singled  out  for  inspec- 
tion.  Thus,  if  ZB  is  dimensioned  (150,150)  and  the  above  values 
are  1,  150,  74,  and  78  respectively,  GRDCHK  will  print  the  depth 
values  of  these  750  elements  of  the  ZB  matrix  and  plot  five 
columnar  profiles  (Fig.  4).   If  this  record  is  omitted,  only 

the  contour  plot  will  be  produced. 
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APPENDIX  E 
ECCOM  OUTPUT  DESCRIPTION  AND  USER  INSTRUCTIONS 


1.  Description  of  Output 

ECCOM  allows  the  user  to  combine  several  selected  ECTRACE  runs  and 
eliminate  specific  rays  as  desired  to  produce  an  improved  horizontal  ray 
path  plot  for  investigation  of  curvature  (Fig.  7  and  Fig.  18).   The  size 
of  the  ECCOM  plot  will  be  the  same  size  as  first  ECTRACE  run  of  the  input 
card  deck  unless  the  user  opts  to  alter  the  scale  factor  (FACT)  to  control 
the  plot  dimensions. 

The  printed  output  is  a  list  of  all  the  rays  processed  for  the  ECCOM 
plot,  giving  each  ray's  initial  elevation  angle  and  heading.   Additionally, 
the  number  of  points  used  to  plot  each  ray  path  is  printed.   This  number 
of  points  corresponds  to  the  number  of  reflections  plus  the  starting  point 
of  the  ray.   Refer  to  Table  VII  for  a  sample  listing  of  data  for  the  ECCOM 
51  plot,  shown  in  Fig.  7. 

2.  User  Instructions 

The  following  JCL  cards  and  parameters  are  required  for  an  ECCOM  run: 

•  JOB  card  -  The  CPU  time  parameter  should  allow  four  seconds  for 
every  55  rays  to  be  plotted. 

•  EXEC  card  -  This  card  must  specify  FORTCLGW  procedure. 

•  Calling  program  -  Following  the 

//FORT. SYS IN  DD  * 
is  the  ECCOM  program  deck.   See  Appendix  L  for  the  program  listing. 

•  Data  cards — Following  the 

//GO. SYS IN  DD  * 
is  the  sequence  of  punched  output  decks,  each  individually  arranged 
as  described  in  Chapter  V,  Section  C.   The  user  may  change  the  title 
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card  (card  #  1)  of  the  first  ECTRACE  punched  deck  to  identify  the 
ECCOM  run  title.   In  addition,  the  user  may  change  the  scale  factor 
(FACT)  on  card  number  4  to  modify  the  plot  size.   FACT  was  calculated 
by  ECTRACE  in  response  to  the  user's  selection  for  the  horizontal 
plot  length  (Chapter  V,  Section  B) .   The  following  are  sample  values 
of  FACT  and  recommended  limits: 

FACT  Horizontal  Plot  Axis  Length  (inches) 

0.333  5  (minimum  recommended) 

1.000  15 

1.267  19  (maximum  recommended) 

If  the  user  desires  a  plot  size  outside  the  range  recommended,  re- 
fer to  Ref.  7. 
Cards  1  through  4  of  each  ECTRACE  deck  after  the  first  are  ignored 
by  ECCOM,  but  the  user  should  verify  that  they  all  reference  the  same 
portion  of  the  same  input  grid-matrix  in  order  for  the  ECCOM  plot  to  be 
meaningful . 
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APPENDIX  F 
ANALYTICAL  VERIFICATION  OF  ECTRACE  HORIZONTAL  RAY  PATHS 


The  accumulation  .of  ray  heading  changes  after  repeated  specular 
reflection  from  a  sloping  sea  bed  describes  a  faceted  horizontal  path 
whose  apparent  curvature  may  in  some  instances  be  analytically  described. 
Harrision  [2]  derived  some  solutions  for  horizontal  ray  paths  over  simple 
bottom  topography  cases.   Most  of  these  cases  may  be  modeled  by  SYNGEN  and 
used  to  test  the  validity  of  the  horizontal  ray  paths  mapped  by  ECTRACE. 
We  chose  an  isovelocity  sound  speed  condition  and  a  conical  seamount  because 
its  surface  geometry  effectively  challenges  the  facet  method  of  topograph- 
ical modeling. 

For  a  conical  seamount  of  arbitrary  slope  but  whose  apex  just  touches 
the  sea  surface,  Harrison  gives  the  ray  paths  in  polar  coordinates  (r,  $ ) 
as 


r  =  r  (1  -  cos  6  cos  <j)  ) 


cosec 


(1  -  cos290cos2(})) 


sin^gCosB, 


2*)2 

(<f>  +  $') 


(F.l) 


where 


initial  distance  from  the  apex 

initial  deviation  angle 

initial  azimuthal  angle 

a  constant  given  by  setting  $  =  0  when  r 


o 


sin4>0cos60 
(1   -   cos2e   cos2(j,  )' 


sin_1(l   -   cos20   cos2cj)   )'5 


(F.2) 
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Figure  F.l.  is  a  diagram  of  a  representative  ray  path,   a  is  the 
deflection  angle  and  represents  the  asymptotic  angle  of  the  shadow  zone  of 
the  ray  fan. 


FIGURE  F.l.   Polar  representation  of  the  horizontal  projection  of  a  ray 
path  for  a  conical  seamount. 


In  the  test  run  (Trial  83,  Fig.  F.2),  the  source  was  positioned  at 
a  range  r  of  50  km  from  the  apex,  bearing  180°,  at  a  depth  just  below 
the  surface.  The  receiver  was  placed  at  the  apex  so  that  for  each  ray 
the  ECTRACE  range  r_  could  be  read  directly  from  the  printed  ray  history 
table  values  for  DIST  TO  REC  and  $  is  the  reciprocal  of  BNG  TO  REC.  A 
zero  gradient  was  used. 
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FIGURE  F.2.  Horizontal  Plot  of  ECTRACE  Trial  for  Harrison  Verification. 
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Two  ray  lines  were  projected  using  6  =  5  degrees  for  ray  fan  1  and 
9  =  10  degrees  for  ray  fan  2.    <j>  ranged  from  270°  for  the  first  ray 

°  Q 

to  90°  for  the  last  ray  in  each  fan,  with  an  increment  of  15°  between  rays. 
Two  rays  were  chosen  for  analysis,  ray  6  (90  =  5  degrees,  (j>0=  -15  degrees) 
and  ray  21  (90  =  10  degrees,  <J>0  =  15  degrees).   The  relevant  portions  of 
their  trace  histories  are  reproduced  in  Table  F.I  and  Table  F.II.   Only 
the  surface  reflection  values  are  used,  where  the  slant  range  DIST  becomes 
the  horizontal  range  r  .   Inserting  9  and  $  values  into  Eq .  (F.l), 

Hi  o  o 

loaded  in  a  programmable  calculator,  we  obtained  a  table  of  values  for  r 
for  both  rays.   The  results  are  shown  in  Table  F.III  and  Table  F.IV  along- 
side the  ECTRACE  values  r  and  the  observed  error  £  =  r  -  r  .   These  same 

E  £• 

results  are  plotted  in  Figures  F.3  and  F.4   (r  vs  $) . 
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TABLE  F.I 


RAY  6  (6  OF  RAY  FAN  1)  TRACE  HISTORY 


REV 

DEPTH 

GRID 

ELEV 

BRNG 

DIST 

NO. 

TYPES 

(KM) 

HEADING 

ANGLE 

TO  REC 

FM  REC 

0 

INIT 

0.001 

345.000 

5.000 

360.000 

50.000 

1 

BOTT 

1.943 

344.500 

-11.825 

11.370 

29.200 

2 

SURF 

0.0 

344.500 

11.825 

22.736 

21.277 

3 

BOTT 

1.153 

343.010 

-16.653 

34.104 

17.329 

4 

SURF 

0.0 

343.010 

16.653 

45.508 

15.170 

5 

BOTT 

0.934 

340.647 

-18.460 

56.912 

14.035 

6 

SURF 

0.0 

340.847 

18.460 

68.425 

13.614 

7 

BOTT 

0.919 

338.231 

-17.062 

79.942 

13.816 

8 

SURF 

0.0 

338.231 

17.062 

91.462 

14.689 

9 

BOTT 

1.096 

336.557 

-12.677 

102.980 

16.464 

10 

SURF 

0.0 

336.557 

12.677 

114.448 

19.713 

11 

BOTT 

1.729 

335.915 

-6.020 

125.912 

25.991 

12 

SURF 

0.0 

335.915 

6.020 

137.458 

40.960 
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TABLE  F.II 
RAY  21  (8  OF  RAY  FAN  2)  TRACE  HISTORY 


REV. 

DEPTH 

GRID 

ELEV 

BRNG 

DIST 

NO. 

TYPES 

(KM) 

HEADING 

ANGLE 

TO  REC 

FM  REC 

0 

INIT 

0.001 

15.000 

10.000 

360.000 

50.000 

1 

BOTT 

2.451 

15.649 

-17.119 

354.386 

36.838 

2 

SURF 

0.0 

15.649 

17.119 

348.770 

29.483 

3 

BOTT 

1.654 

17.253 

-23.468 

343.153 

24.866 

4 

SURF 

0.0 

17.253 

23.468 

337.518 

21.761 

5 

BOTT 

1.303 

20.Q43 

-28.702 

331.883 

19.590 

6 

SURF 

0.0 

20.043 

28.702 

326.243 

18.046 

7 

BOTT 

1.128 

24.047 

-32.594 

320.603 

16.949 

8 

SURF 

0.0 

24.047 

32.594 

315.014 

16.200 

9 

BOTT 

1.046 

28.867 

-33.973 

309.425 

15.723 

10 

SURF 

0.0 

28.867 

33.973 

303.767 

15.480 

11 

BOTT 

1.029 

34.081 

-33.685 

298.105 

15.459 

12 

SURF 

0.0 

34.081 

33.685 

293.479 

15.661 

13 

BOTT 

1.071 

38.566 

-31.042 

286.852 

16.097 

14 

SURF 

0.0 

38.566 

31.042 

281.203 

16.802 

15 

BOTT 

1.187 

41.986 

-26.620 

275.55 

17.834 

16 

SURF 

0.0 

41.986 

26.620 

269.889 

19.295 

17 

BOTT 

1.420 

44.169 

-20.799 

264.224 

21.345 

18 

SURF 

0.0 

44.169 

29.799 

258.538 

24.279 

19 

BOTT 

1.904 

45.278 

-14.066 

252.851 

28.620 

20 

SURF 

0.0 

45.278 

14.066 

247.160 

35.467 

21 

BOTT 

3.160 

45.674 

-6.733 

241.469 

47.511 

22 

SURF 

0.0 

45.674 

6.733 

235.781 

73.529 
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TABLE  F 

.III 

Trial 

83 

Ray  6 

$ 

r 

rE 

e 

0.0° 

50.0  km 

50.0  km 

0.0km 

-22.736 

21.262 

21.277 

-0.015 

-45.408 

15.176 

15.170 

0.006 

-68.425 

13.616 

13.614 

0.002 

-91.462 

14.713 

14.689 

0.024 

114.49 

19.822 

19.713 

0.109 

137.46 

41.571 

40.960 

0.611 

Ray   7        r^  4.339   km  ($  =   d? .   constant) 

E 
ce 


-   r^        =0.019   km 
ce  E 

ce 
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TABLE  F 

.IV 

Trial  83 

Ray  21 

$ 

r 

rE 

£ 

0.0° 

50.0  km. 

50.0  km 

0.0  km 

11.230 

29.470 

29.483 

-0.013 

22.482 

21.745 

21.761 

-0.016 

33.757 

18.026 

18.046 

-0.020 

44.986 

16.178 

16.200 

-0.022 

56.233 

15.548 

15.480 

-0.022 

67.521 

15.643 

15.661 

-0.018 

78.797 

16.791 

16.802 

-0.011 

90.111 

19.309 

19.295 

0.014 

101.46 

24.365 

24.279 

0.086 

112.84 

35.822 

35.467 

0.355 

124.219 

75.857 

73.529 

2.238 

Ray   20   r_        =   8.687   km  ($  =   0,    constant) 

E 
ce 

r        -  r„        =  -0.005  km 
ce  E 

ce 
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FIGURE  F.3.   Harrison 


verificaton  plot  of  range  vs  $  for  ray  6. 
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FIGURE  F.4.   Harrison  Verification  Plot  of  range  vs  $  for  ray  21 
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These  results  show  very  close  agreement  up  to  a  distance  well  past 
the  closest  point  of  approach,  and  only  begin  to  diverge  slightly  when  the 
rays  approach  their  asymptotic  values  for  $.   Investigation  of  other  rays 
of  the  same  trial  supports  this  observation.   Clearly  the  fact  that  the 
plane  facets  do  not  all  have  the  same  slope  as  the  cone  they  collectively 
represent,  but  vary  as  much  as  4%  between  their  extremes,  explains  most 
of  the  error.   The  validity  of  the  model  decreases  near  the  apex,  where 
fewer  data  points  and  hence  plane  facets  are  used  to  approximate  the 
conical  surface.   In  addition,  the  progressively  shallower  grazing  angles 
assumed  by  the  rays  after  CPA  cause  the  horizontal  paths  to  become  more 
sensitive  to  this  facet  effect. 

Reference  13  also  provides  an  equation  for  the  closest  approach  of 
the  ray  fan  envelope  to  the  seamount  center,  developed  by  letting 
4>0  =  0  in  Eq.  (F.l)  .   Thus 


r   =  r   sin  t) 
ce    o      q 


=  8.682  km.  (F.3) 


on 


Ray  21  of  Trial  83  was  the  radial  ray  (4>0  =  0°)  and  its  CPA 
the  printout  was  8.687  km,  a  difference  of  only  four  meters  after  11 
reflections  from  the  bottom  model.   The  ray  should  have  reversed  its 
heading  exactly  180  degrees  after  CPA,  but  actually  experienced  a  174 
degree  heading  change,  revealing  the  sensitivity  of  the  model  to  radially- 
directed  ray  paths. 

Lastly,  Ref.  13  defines  the  asymptotic  angle  of  the  shadow  zone  of 
the  ray  fan  (Fig.  F.l).   A  plot  of  the  final  azimuthal  angle  $     versus 
6Q  for  ray  fan  2  (Fig.  F.5)  of  the  test  run  exhibits  an  interpolated 
value  for  a  of  40  degrees  for  <|>  =  29  degrees.   Testing  this  value  with 

95 


90_ 


a  = 


Asymptotic  deflection 
angle  for  ray  fan  of 
S   =10  degrees 


(j)f  VS  <j) 


a  - 


-45 


i 
45 


-  a 


--90 
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the  formula 

o  1/ 

a  =  4<Xsin29n)  /3 

=  39.187  degrees  (F.4) 

shows  close  agreement  with  the  theoretical  value.   Figure  F.6,  an  ECCOM 
composite  for  the  seamount,  is  a  graphic  example  of  the  dependence  of 
a  on  the  initial  elevation  angle. 

It  is  worth  noting  that  the  analytical  expressions  for  horizontal 
projections  of  ray  paths  apply  only  to  exact  functional  bottom  topographies 
and  that  the  curved  three  dimensional  ray  path  surfaces  only  approximate 
an  actual  surface  of  vertical  planes  joined  at  bottom  reflection  points. 
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ECC8H    28      »««    SYNTHETIC    SEflM8UNT 
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FIGURE  F.6.  ECCOM  Trial  28  showing  the  effect  of  a  synthetic  conical 
seamount  on  ray  propagation.   The  ray's  initial  headings 
varied  from  20°  to  58°  and  the  initial  elevation  varied 
from  4°  to  18°.   Selected  rays  from  five  ECTRACE  runs 
(100  rays)  were  combined. 
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APPENDIX  G 
TABLES 
TABLE  I 
A  SUMMARY  OF  THE  PROGRAMS  AND  THEIR  COMPONENT  SUBPROGRAMS 


I.   GENBOT  -  plots  source  contour  data  and  generates  depth  matrix 

A.  CORNER  -  calculates  boundary  coordinates 

B.  GEOPLT  -  plots  geometric  grid  lines 

C.  GEODST  -  solves  geodetic  triangle 

D.  RAIN  1/2/3  -  perform  sorting  and  interpolation  (from  SSP3) 

II.   SYNGEN  -  produces  synthetic  sea  bed  depth  matrices 

III.   G3DP  -  produces  perspective  surface  plot  of  depth  matrix 

A.  GRDSCT  -  loads  working  surface  plot  of  depth  matrix 

B.  PLT3D1  -  draws  surface  plot  (from  NPS) 

IV.   GRDCHK  -  allows  inspection  of  depth  matrix  and  draws  contour  plot 
of  same 

A.  GRDSCT  -  see  above 

B.  CORNER  -  see  above 

C.  GEOPLR  -  same  as  GEOPLT  but  for  a  rotated  plot 

D.  GEODST  -  see  above 

E.  CONTUR  -  draws  a  contour  plot  of  depth  matrix  (from  NPS) 

V.   ECTRACE  -  a  collection  of  subprograms  which  perform  ray  tracing 

A.  TRACER  -  controls  all  subprograms  and  manages  ray  traces 

B.  GRDSCT  -  see  above 

C.  IDTSUB  -  identifies  subscripts  of  depth  matrix 

D.  DEEP  -  calculates  depth  at  given  coordinates 

E.  IDPROF  -  identifies  water  mass 

F.  NUPROF  -  calculates  water  mass  parameters 

G.  CHNLIM  -  calculates  sound  channel  limits 

H.  CONTAC  -  calculates  point  of  bottom  contact 

I.  BOANG  -  calculates  direction  and  grazing  angle  of  bottom 

reflected  ray 

J.  FBLOSS  -  calculates  bottom  loss 

K.  ANGPRT  -  prints  ray  fan  bottom  bounce  table 

L.  BNGPLT  -  initializes  plotting  and  draws  plot  borders 

M.  BDTPLT  -  draws  bathymetry  profile 

N.  RNGPLT  -  draws  vertical  plot 

0.  T2DPLT  -  draws  horizontal  plot 

P.  SSPPLT  -  draws  sound  speed  profiles 

Q.  ENDPLT  -  draws  plot  legend  and  terminates  plotting 

(L  through  Q  call  subroutines  from  the  Versatec  software  package.) 
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TABLE  II 
INPUT  PARAMETERS  TABLEAU  (COMPUTER  PRINT  OUTPUT) 


ECTRACE  TRIAL  67  ***  SYNTHETIC  BASIN  *** 


SOUND  SPEED  PROFILES 


ZC 


C(ZC) 


GRAD 


0.0  1461.08 

50.00  1461.00 

60.00  1447.00 

70.00  1445.00 

90.00  1447.00 

500.00  1454.30 

1050.00  1461.00 

1500.00  1470.00 

4000.00  1511.70 

NWMR  (Northwest  Mohns  Ridge) 


0.0 

■1.40000 
-0.20000 
0.10000 
0.01780 
0.01218 
0.02000 
0.01668 
0.01668 


ZC 


C(ZC) 


GRAD 


0.0  1490.00 

50.00  1482.00 

200.00  1474.50 

375.00  1474.00 

750.00  1466.00 

975.00  1465.00 

1250.00  1467.00 

1500.00  1470.00 

4000.00  1511.70 

SEMR  (Southeast  Mohns  Ridge) 


■0.16000 

■0.05000 

■0.00286 

-0.02133 

-0.00444 

0.00727 

0.01200 

0.01668 

0.01668 


ZC 


C(ZC) 


GRAD 


0.0 
1500.00 
4000.00 


1470.00 
1470.00 
1511.70 


0.0 

0.01668 

0.01668 


ISOV  (Isovelocity) 


(Table  II  continued) 
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continued  TABLE  II 


FRONT /BOUNDARY  INTERSECTIONS 

FRONT  FRONT 

ENTRY  EXIT  ENTRY  EXIT 

(  0.0,  50.0)   (100.0,150.0)         (  50.0,  0.0)      (150.0,  100.0) 

GRID  FILE  CALLED   ***************************************************** 

*  SYNTHETIC  SEA  BED  * 

*  COMPUTER  GENERATED   BASIN  * 
***************************************************** 

HAD  BEEN  READ.   CENTER  LAT/LONG  IS   70.07059/  15.85032  AT 

(75.000,  75.000)  KM.  RANGES  ARE   0.0  TO  150.0  KM  IN  Y  (22801  ITEMS). 

THE  FOLLOWING  ADDITIONAL  PARAMETERS  WERE  INPUT: 

SOURCE  DATA:   FREQUENCY  150.00  HZ, GRID  COORDINATES  (  40.00,40.00), 
DEPTH   0.850  KM.   8  RAY  FANS  WILL  BE  TRACED,  FROM  -30.0  TO  30.0  DEGREES 
ELEVATION,  INCREMENTING  EVERY   8.6  DEG.   1  RAYS  MAKE  A  RAY  FAN  FROM   10.0 
TO   10.0  DEGREES  TRUE  HEADING,  INCREMENTING  EVERY   0.0  DEGREES. 

RECEIVER  DATA:   GRID  COORDINATES  (125.0,   40.0),  DEPTH  1.749  KM. 

GRID  ARRAY:   MAXIMUM  ARRAY  SIZE  (151 ,151) .DEPTH  VALUES  RESOLVED  BY   1.0 
KM.  SPACING. 

TRACE  PARAMETERS:   PROPAGATION  WILL  CONTINUE  UNTIL  A  BOTTOM  BOUNCE  LOSS 

OF  150.0  DB  OR   21  BOTTOM  CONTACTS  EXCEEDED  FOR  EACH 
RAY.  PATH  INCREMENT  IS   100.0  M.  POSITION  ERROR  AFTER 
BOTTOM  CONTACT  WITHIN   0.05000  KM.   CHANNEL  RANGE 
LIMIT  200.0  KM.  SMALL  ANGLE  2.0  DEG. 


PLOT  OPTION  SELECTED. 

PLOT  PARAMETERS: 

X 

Y 

Z 

R 

AXES: 

0.0 

0.0 

-1.0 

-10.0 

10.0 

10.0 

2.0 

7.00 

150.0 

150.0 

6.0 

123.00 

SCALING:   LENGTH  4  IN.,  FACTOR=0. 2667 
EACH  SEGMENT  REPRESENTS  10  POSITION  CALCULATIONS 


(program  parameters  are  user  supplied  via  DATA  cards) 
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TABLE  III 

Ray  History  Tableau  ECTRACE  Trail  67  (computer  print  output) 
Ray  8   (of  Ray  Fan  8)  TRACE  HISTORY 


REV.  NO 

.   TYPE 

X(KM) 

Y(KM) 

DEPTH (ft 

*)             HEADING 

ELEV.  ANG. 

0 

INIT 

40.000 

40.000 

0.850 

10.000 

30.000 

1 

BOTT 

40.722 

44.097 

3.161 

12.533 

-20.323 

2 

SURF 

42.397 

51.631 

0.0 

12.533 

21.104 

3 

BOTT 

44.611 

61.591 

4.043 

14.157 

-14.658 

4 

SURF 

47.628 

73.549 

0.0 

14.157 

17.615 

5 

BOTT 

50.933 

86.652 

4.344 

15.556 

-14.938 

6 

SURF 

54.384 

99.050 

0.0 

15.556 

18.425 

7 

BOTT 

57.154 

109.001 

3.628 

17.472 

-22.313 

8 

SURF 

59.508 

116.479 

0.0 

17.472 

26.133 

9 

BOTT 

61.260 

122.046 

2.780 

20.496 

-31.868 

10 

SURF 

62.747 

126.024 

0.0 

20.496 

33.698 

11 

BOTT 

63.904 

129.119 

2.190 

25.415 

-41.798 

12 

SURF 

64.930 

131.279 

0.0 

25.415 

42.680 

13 

BOTT 

65.778 

133.064 

1.823 

33.460 

-51.533 

14 

SURF 

66.569 

134.260 

0.0 

33.460 

51.978 

15 

BOTT 

67.260 

135.305 

1.606 

46.953 

-60.129 

16 

SURF 

67.933 

135.934 

0.0 

46.953 

60.372 

17 

BOTT 

68.551 

136.511 

1.490 

69.313 

-65.773 

18 

SURF 

69.182 

136.749 

0.0 

69.313 

65.930 

19 

BOTT 

69.787 

136.977 

1.452 

97.778 

-66.103 

20 

SURF 

70.428 

136.890 

0.0 

97.778 

66.244 

21 

BOTT 

71.073 

136.802 

1.482 

121.284 

-61.160 

22 

SURF 

71.771 

136.378 

0.0 

121.284 

61.347 

23 

BOTT 

72.510 

135.929 

1.587 

135.723 

-53.031 

24 

SURF 

73.338 

135.080 

0.0 

135.723 

53.338 

25 

BOTT 

74.265 

134.129 

1.786 

144.372 

-43.592 

26 

SURF 

75.339 

132.631 

0.0 

144.372 

44.152 

27 

BOTT 

76.618 

130.845 

2.127 

149.526 

-34.065 

28 

SURF 

78.158 

128.228 

0.0 

149.526 

35.164 

29 

BOTT 

80.046 

125.020 

2.664 

153.024 

-26.290 

30 

SURF 

82.376 

120.443 

0.0 

153.024 

26.244 

31 

BOTT 

85.445 

114.414 

3.452 

155.238 

-19.001 

32 

SURF 

89.145 

106.392 

0.0 

155.238 

20.346 

33 

BOTT 

93.855 

96.179 

4.247 

156.763 

-15.037 

34 

SURF 

98.821 

84.615 

0.0 

156.763 

18.321 

35 

BOTT 

103.807 

73.002 

4.261 

158.297 

-17.529 

36 

SURF 

107.968 

62.549 

0.0 

158.297 

20.418 

37 

BOTT 

111.257 

54.284 

3.486 

160.557 

-25.605 

38 

SURF 

113.497 

47.940 

0.0 

160.557 

28.126 

39 

BOTT 

115.211 

43.084 

2.713 

164.103 

-35.243 

40 

SURF 

116.231 

39.501 

0.0 

164.103 

36.331 

41 

BOTT 

117.049 

36.630 

2.188 

169.812 

-45.355 

BOTTOM  LOSS  GATE 

EXCEED.   TERMINATINC 

;  RAYTRACE  AFTER 

3311 

SEGMENTS 

.   PLOT  BASED  ON  367 

POINTS, 

41  REVERSALS  ANI 

)  21 

BOTTOM  BOUNCES. 

MINIMUM  DISTANCE  TO 

RECEIVER,   8.57 

KM, 

OCCURS  AT  TIME   166.244 

(TABLE  III  continued) 
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contir 

iued  TABLE 

III 

BRING 

DIST 

WTR 

DEPTH .. 

LIMITS 

TO  SRC* 

FM  SRC* 

TIME 

MASS 

LYR 

MIN 

MAX 

0.0 

0.0 

0.0 

SEMR 

5 

0.0 

•kit-kick 

190.00 

4.759 

3.270 

SEMR 

8 

0.0 

***** 

191.646 

11.906 

8.802 

SEMR 

1 

0.0 

***** 

192.056 

23.308 

16.265 

SEMR 

9 

0.0 

7.093 

192.809 

34.416 

24.931 

SEMR 

1 

0.0 

7.093 

193.189 

48.043 

34.520 

SEMR 

9 

0.0 

7.526 

193.690 

60.783 

43.604 

SEMR 

1 

0.0 

7.526 

193.961 

71.155 

51.036 

NWMR 

8 

0.0 

***** 

194.310 

78.932 

56.802 

NWMR 

1 

0.0 

***** 

194.572 

84.778 

61.279 

NWMR 

8 

0.0 

***** 

194.812 

88.984 

64.648 

NWMR 

1 

0.0 

***** 

195.015 

92.279 

67.440 

NWMR 

8 

0.0 

***** 

195.276 

94.626 

69.600 

NWMR 

1 

0.0 

***** 

195.483 

96.573 

71.568 

NWMR 

8 

0.0 

***** 

195.741 

97.936 

73.103 

NWMR 

1 

0.0 

***** 

195.962 

99.130 

74.624 

NWMR 

8 

0.0 

***** 

196.234 

99.921 

75.814 

NWMR 

1 

0.0 

***** 

196.480 

100.647 

77.058 

NWMR 

7 

0.0 

***** 

196.784 

101.058 

78.109 

NWMR 

1 

0.0 

***** 

197.074 

101.451 

79.258 

NWMR 

7 

0.0 

***** 

197.435 

101.559 

80.307 

NWMR 

1 

0.0 

***** 

197.796 

101.669 

81.494 

NWMR 

7 

0.0 

***** 

198.245 

101.483 

82.550 

NWMR 

1 

0.0 

***** 

198.721 

101.290 

83.914 

NWMR 

8 

0.0 

***** 

199.322 

100.758 

85.183 

NWMR 

1 

0.0 

***** 

200.003 

100.176 

86.827 

NWMR 

8 

0.0 

***** 

200.882 

99.146 

88.467 

NWMR 

1 

0.0 

***** 

201.954 

97.956 

90.673 

NWMR 

8 

0.0 

***** 

203.388 

96.130 

93.144 

NWMR 

1 

0.0 

***** 

205.222 

93.996 

96.346 

SEMR 

8 

0.0 

***** 

207.779 

90.926 

100.139 

SEMR 

1 

0.0 

***** 

211.412 

87.232 

105.373 

SEMR 

8 

0.0 

***** 

216.510 

82.606 

111.713 

SEMR 

1 

0.0 

***** 

223.790 

77.898 

110.854 

SEMR 

9 

0.0 

7.469 

232.820 

73.832 

128.722 

SEMR 

1 

0.0 

7.469 

242.651 

71.917 

137.727 

SEMR 

9 

0.0 

8.687 

251.646 

71.615 

145.718 

SEMR 

1 

0.0 

8.687 

258.664 

72.722 

152.236 

ISOV 

2 

0.0 

***** 

263.834 

73.929 

157.358 

ISOV 

1 

0.0 

***** 

267.652 

75.297 

161.318 

ISOV 

2 

0.0 

***** 

270.375 

76.237 

164.413 

ISOV 

1 

0.0 

***** 

272.504 

77.134 

166.955 

ISOV 

2 

0.0 

***** 

*For  this  run  range  and  bearing  of  the  rays  were  referenced  to  the 
source.   As  explained  in  Chapter  V,  in  the  current  version  of 
ECTRACE,  these  values  are  now  referenced  from  the  receiver  (REC) . 

(TABLE  III  continued) 
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continued  TABLE   III 


c 

GRAZ 

BTM 

BTM 

(KM/S) 

ANGLE 

ANGLE 

LOSS 

LPS 

1.46556 

1.49770 

24.081 

4.678 

4.1 

0 

1.49000 

1.51242 

16.731 

3.392 

6.1 

0 

1.49000 

1.51743 

14.465 

2.756 

7.2 

0 

1.49000 

1.50549 

19.479 

3.935 

10.0 

1 

1.46100 

1.49135 

27.785 

5.015 

15.1 

1 

1.46100 

1.48151 

37.224 

5.651 

22.2 

1 

1.46100 

1.47539 

46.949 

6.015 

30.0 

0 

1.46100 

1.47177 

56.159 

6.182 

40.9 

0 

1.46100 

1.46980 

63.522 

6.272 

50.9 

0 

1.46100 

1.46904 

66.735 

6.262 

60.0 

0 

1.46100 

1.46963 

64.227 

6.248 

70.9 

0 

1.46100 

1.47145 

57.381 

6.141 

80.9 

0 

1.46100 

1.47477 

48.418 

6.040 

89.9 

0 

1.46100 

1.48045 

38.786 

5.639 

97.3 

0 

1.46100 

1.48941 

30.571 

5.162 

103.0 

2 

1.49000 

1.50255 

22.176 

4.155 

106.6 

0 

1.49000 

1.51582 

16.287 

2.886 

108.3 

0 

1.49000 

1.51606 

16.298 

2.906 

110.1 

1 

1.49000 

1.50313 

22.375 

4.256 

113.8 

0 

1.47000 

1.49023 

30.985 

5.212 

119.5 

0 

1.47000 

1.48147 

40.627 

5.821 

127.3 

0 
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TABLE  IV 

Ray  fan  Summary  Tableau  for  ECTRACE  Trail  67  ray  number  8.  (computer  print 
output) 

ANGLES  BEFORE  BOUNCE.   INITIAL  ELEVATION  ANGLE  =  30.00 


HEADINGS : 

10.00 

BOUNCE 

1 

27.78 

2 

18.77 

3 

13.96 

4 

16.60 

5 

23.63 

6 

32.51 

7 

42.10 

8 

51.66 

9 

60.18 

10 

65.80 

11 

66.12 

12 

61.15 

13 

52.96 

14 

43.36 

15 

34.76 

16 

25.30 

17 

17.50 

18 

15.03 

19 

19.09 

20 

26.64 

21 

35.75 

AT  BOUNCE. 

INITIAL  ELEVATION  ANGLE 

HEADINGS : 

10.00 

BOUNCE 

1 

3160.6 

2 

4042.9 

3 

4343.8 

4 

3627.9 

5 

2780.2 

6 

2190.1 

7 

1823.4 

8 

1606.3 

9 

1490.1 

10 

1452.0 

11 

1481.7 

12 

1586.9 

13 

1785.7 

14 

2126.5 

15 

2663.6 

16 

3451.5 

17 

4246.8 

18 

4261.3 

19 

3486.1 

20 

2713.1 

21 

2187.8 

=  30.00 
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TABLE  VI 

JAY  FAN  SYMBOLS  -  used  to  denote  location  of  bottom  reflections  in  the 
£CTRACE  plots.   Thirteen  different  symbols  are  used  before  repetition. 
U.1  rays  of  the  same  ray  fan  (initial  elevation  angle)  will  indicate  a 
jottom  reflection  with  the  same  symbol.   Each  ECTRACE  plot  lists  the 
symbols  used  for  the  horizontal  and  vertical  plots. 


INITIAL 

ELEVATION 

ANGLES 

RnTFRN 

SYMBOL 

THETA 

1 

□ 

-10.000 

2 

O 

-9.000 

3 

A 

-8.000 

4 

-7.000 

5 

X 

-6.000 

6 

<!> 

-5.000 

7 

+ 

-4.000 

8 

X 

-3.000 

9 

z 

-2.000 

10 

Y 

-1.000 

11 

X 

0.000 

12 

* 

1.000 

13 

X 

2.000 

14 

□ 

3.000 

15 

0 

4.000 

16 

A 

5.000 

17 

+ 

6.000 

18 

X 

7.000 

19 

o 

8.000 

20 

+ 

9.000 

21 

X 

10.000 

HEADINGS  90.000  +/-  0.000 
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TABLE  VII 

ECCOM  printed  output.   The  rays  which  are  plotted  are  listed,  and  their 
initial  heading,  initial  elevation,  and  the  number  of  points  (initial  plus 
reflections)  are  indicated.   In  this  run  the  plot  produced  is  shown  in 
Fig.  7. 


ECCOM  51  ***  SYNTHETIC   BASIN  *** 
SYNTHETIC  SEA  BED 
COMPUTER  GENERATED 


BASIN 

XDV        YST 

YDV       CNX 

CNY 

CNLAT     CNLON 

10.000      0.0 

10.000    75.000 

75.000 

0.0      0.0 

INTIAL 

INTIAL 

ELEVATION 

HEADING 

POINTS 

SCALE  FACT0R=  0.3 

40.000 

40.000 

21 

40.000 

44.000 

21 

40.000 

48.000 

21 

38.000 

42.000 

21 

38.000 

46.000 

21 

36.000 

40.000 

21 

36.000 

44.000 

21 

36.000 

48.000 

21 

34.000 

42.000 

21 

34.000 

46.000 

21 

60.000 

40.000 

21 

60.000 

44.000 

21 

60.000 

48.000 

21 

51.333 

42.000 

21 

51.333 

46.000 

21 

42.667 

40.000 

21 

42.667 

44.000 

21 

42.667 

48.000 

21 

34.000 

42.000 

21 

34.000 

46.000 

21 

60.000 

5.000 

21 

60.000 

9.000 

21 

60.000 

13.000 

21 

62.000 

7.000 

21 

62.000 

11.000 

21 

64.000 

5.000 

21 

64.000 

9.000 

21 

64.000 

13.000 

21 

66.000 

7.000 

21 

66.000 

11.000 

21 
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APPENDIX  H 
FIGURES 


FIGURE  1.  Program  and  Data  Utilization  Flowchart 

FIGURE  2.  Contour  Plot  of  Lofoten  Basin  (GENBOT) 

FIGURE  3.  Horizontal  Ray  Plot  Over  Mohns  Ridge  (ECTRACE) 

FIGURE  4.  Contour  Plot  of  Mohns  Ridge  (GRDCHK) 

FIGURE  5.  Complete  ECTRACE  Plot  for  a  Synthetic  Basin 

FIGURE  6.  Contour  Plot  of  a  Synthetic  Basin  (GRDCHK) 

FIGURE  7.  ECCOM  Plot  for  a  Synthetic  Basin 

FIGURE  8.  NRL  N.  Atlantic  Contour  Chart 

FIGURE  9.  Perspective  Surface  Plot  -  Seamount  (SYNGEN) 

FIGURE  10.  Perspective  Surface  Plot  -  Basin  (SYNGEN) 

FIGURE  11.  Perspective  Surface  Plot  -  Mohns  Ridge  (SYNGEN) 

FIGURE  12.  Contour  Plot  of  a  Synthetic  Seamount  (GRDCHK) 

FIGURE  13.  Complete  ECTRACE  plot  for  a  Synthetic  Seamount. 

FIGURE  14.  Horizontal  Ray  Plot  over  a  Synthetic  Seamount  (ECTRACE) 

FIGURE  15.  Vertical  Ray  Plot  in  a  Basin  (ECTRACE) 

FIGURE  16.  Vertical  Ray  Plot  over  a  Seamount  (ECTRACE) 

FIGURE  17.  Horizontal  Ray  Plot  in  a  Basin  (ECTRACE) 

FIGURE  18.  ECCOM  Plot  for  a  Synthetic  Wedge 
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FIGURE  1.   Program  and  Data  Utilization  Flowchart 
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FIGURE  2.   A  contour  plot  of  the  eastern  Lofoten  Basin  produced  by 
GENBOT  from  magnetic  tape  data.   The  rectantular  labeled 
axis  denotes  the  region  for  which  GENBOT  generated  a 
matrix  representing  evenly  spaced  depth  values.   A 
geographic  coordinate  grid  is  superimposed. 
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FIGURE  3.   ECTRACE  TRIAL  82C  performed  over  a  depth-grid  matrix  of 
the  Mohns  Ridge.   In  this  run  three  ray  fans  are  used 
with  eight  rays  per  ray  fan.   This  may  be  overlaid  on 
a  contour  plot  of  this  area  (Fig.  2). 
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A  GRDCHK  plot  made  on  a  matrix  generated  from  contour  data. 
In  this  example,  rows  1  through  150  of  columns  74  through  78 
were  selected  for  profile  analysis  and  printing. 
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FIGURE  5.   ECTRACE  TRIAL  67,  performed  in  a  synthetic  parabolic  basin 

(FIGURE  6) .   In  this  run  eight  ray  fans  are  used  with  initial 
elevation  angles  varying  between  -30  and  +30  degrees.  Since 
there  is  only  one  ray  in  each  ray  fan,  the  initial  heading 
for  the  rays  will  be  10  degrees.   The  vertical  plot  reveals 
CZ,  sofar  ducted,  and  SRBR  rays. 
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FIGURE  6.   A  graphical  analysis  performed  by  GRDCHK  of  a  synthetic 

parabolic  basin  (Fig.  10)  created  by  SYNGEN.  The  upper 

plot  is  a  series  of  displaced  depth  curves  of  selected 
columns  of  the  matrix. 
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FIGURE  7.   Convergence  in  a  parabolic  basin  (Fig.  10).   This  is  an 
ECCOM  product  using  punched  output  from  three  separate 
ECTRACE  runs.   TABLE  VII  identifies  the  ray  headings, 
elevation  angles,  and  number  of  points. 
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FIGURE  8.   A  reduced  representation  of  a  chart  titled  "Bathymetry  of 

the  Norwegian-Greenland  and  Western  Barents  Seas"  published 
by  NRL  (Washington,  D.C.).   The  actual  chart  measures  4  by 
3  feet  depicting  finely  grained  contours.   GENBOT  uses 
bathymetry  data  files  of  this  area  to  generate  the  matrix 
grids  used  by  ECTRACE. 
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SYNTHETIC    SER    BED 
COMPUTER    GENERATED 


SERMGUNT 


FIGURE  9.   A  perspective  surface  plot  produced  by  G3DP  of  a  conical 
seamount  generated  by  SYNGEN. 
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SYNTHETIC  SER  BED 
COMPUTER    GENERATED 


BASIN 


FIGURE   10. 


A  perspective  surface  plot  of  a  synthetic  parabolic  basin 
produced  by  SYNGEN. 
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MOHNS  RIDGE 
CENTER  LflT/LGN  7 1 -37N/000-02W 

FIGURE  11.   A  perspective  surface  plot  produced  by  G3DP  of  a  depth 
matrix  representing  a  portion  of  the  Mohns  Ridge,  see 
Fig.  2. 
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FIGURE  12.   A  graphical  analysis  performed  by  GRDCHK  of  a  depth  matrix 
representing  a  conical  seamount.   The  upper  portion  is  a 
series  of  displaced  plots  of  selected  columns  of  the  matrix, 
The  lower  portion  is  a  contour  plot.   The  seamount  was 
created  by  SYNGEN,  see  Fig.  9. 
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FIGURE  13.   An  ECTRACE  plot  of  two  ray  fans  in  a  synthetic  parabolic 
basin  (FIGURE  10)  generated  by  SYNGEN.   The  ray  fans  have 
initial  elevation  angle  of  30  and  45  degrees  (positive 
angle  downward) .   Each  ray  fan  contains  two  rays  with 
initial  headings  360  and  020  degrees,  respectively.  The 
plot  was  terminated  after  each  ray  had  accumulated  23 
bottom  bounces. 
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40. CO  50.00 

FIGURE   17. 


BO. 00  70.00  80.00  90.00  100.00         110.00         150.00         132. CO 

An  ECTRACE  horizontal  plot  product  of  a  single  wide  ray  fan 
in  a  synthetic  basin.   All  rays  share  the  same  bottom  con- 
tact symbol.   The  number  at  the  last  symbol  identifies  the 
ray's  order  in  the  printed  and  punched  history. 
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FIGURE  18.   A  horizontal  plot  produced  by  ECCOM  of  four  ray  fans  propa- 
gating over  a  six  degree  wedge.   The  punched  output  of  three 
ECTRACE  runs  was  combined  and  then  every  other  ray  data  his- 
tory card  sequence  was  removed.   There  are  a  total  of  40  rays 
depicted  varying  in  initial  elevation  angles  from  20  to  26 
degrees  in  two  degree  increments.   The  initial  headings  vary 
from  31  to  49  in  two  degree  increments.   The  main  purpose  of 
ECCOM  is  to  display  the  horizontal  curvature  effect  of  a 
large  number  of  rays. 
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APPENDIX  0 
FORTCLGW  JCL 


The  job  control  language  (JCL)  card  for  FORTCLGW  used  in  several  of 
the  programs  was 

//EXEC  FORTCLGW.  REG  I  ON.  GO  2  5  OK   . 

This  JCL  card  identifies  and  calls  the  following  JCL: 

XXFORTCLGW  PROC  DEST=A, TMSL=SF 

XXFORT   EXEC   PGM=IEYFORT,REGION=180K 

XXSYSPRINT  DD  SYSOUT=&DEST,DCB= (RECFM=FBA^LRECL=120,BLKSIZE=3360) 

SSSYSLIN  DD  UNIT=SYSDA,SPACE-(CYL, (1 , 1) ) ,DISP=( ,PASS) ,DSN=&&SYSLIN, 

XX  DOE= (RECFM=  FB , LRECL=80 , BKJSUZE=800) 

XXLINK  EXEC   PGM=IEWL,REGI0N=180K,PARM= 'MAP ,LIST ' ,COND= (4 ,LT) 

XXSYSLIB  DD   DISP=SHR,DSN=SYS1.F0RTLIB 

XX       DD   DISP=SHR,DSN=SYS1.MPSLIB 

XX       DD  DISP=SHR,DSN=SYS3.IMSL.&IMSL 

XX       DD  DISP=SHR,DSN=SYS1.VTECPL0T 

XXSYSLMOD  DD  DSN=&T (PROGRAM) ,UNIT=SYSDA,DISP=( , PASS) , 

XX   SPACE=(CYL, (3,1,1)) 

XXSYSLIN   DD   DSN=*. FORT. SYSLIN,DISP= (OLD, DELETE) 

XX       DD   DISP=SHR,DSN=SYS1.PR0CLIE(VMAPP) 

XX        DD  DDNAME=SYSIN 

XXSYSPRINT  DD    SYSOUT=&DEST,DOE=(RECFM=FBA,LRECL=121 ,BLKSIZE=1210) , 

XX   SPACE=(CYCL,(l,l)),UNIT=(SYSOUT,SEP=SYSLIN) 

XXSYSUT1   DD   UNIT=(SYSDA,SEP=(SYSLIN,SYSLMOD,SYSLIB)), 

XX   SPACE=(TRK,(10,5)) 

XXPLOT  EXEC  PGM=IEVMAPP,C0ND=(4LT) 

XXSTEPLIB  DD   DSN=SYS1 . VTECPLOT,DISP=SHR 

XXPLOTLOG   DD   SYSOUT=A 

XXVECTR1   DD  DISP= (OLD, DELETE) ,DSN=&&VECTRl 

XXSYSVECTR  DD   DISP= (OLD , DELETE) ,DSN=&&VECTR2 

XXSYSVECTOR  DD  SYSOUT= (A, , 5555) 

XXVECTTAPE   DD  DUMMY 
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